https://femiguez.github.io/nlraa-docs/nlraa-Oddi-LFMC.html

https://onlinelibrary.wiley.com/doi/pdfdirect/10.1002/ece3.5543

# install.packages("emmeans")

library(nlme)
library(nlraa)
library(ggplot2)
library(emmeans) # emmeans::contrast(), emmeans::emmeans()
library(tidyverse)

Steps in fitting a nonlinear mixed model

Part I.

data(lfmc)
glimpse(lfmc)
Rows: 247
Columns: 6
$ leaf.type <fct> M. spinosum, M. spinosum, M. spinosum, S. bracteolactus, S. bracteolactus, S. bracteolactus, Gras…
$ time      <int> 1, 1, 1, 1, 1, 1, 1, 1, 13, 13, 13, 13, 13, 13, 13, 13, 13, 26, 26, 26, 26, 26, 26, 26, 26, 26, 3…
$ plot      <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
$ site      <fct> P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P, P…
$ lfmc      <dbl> 238.4, 269.8, 325.2, 326.7, 257.3, 279.3, 91.9, 95.0, 215.8, 216.7, 205.3, 208.5, 225.7, 246.8, 8…
$ group     <fct> SM plot 1, SM plot 1, SM plot 1, SS plot 1, SS plot 1, SS plot 1, GE plot 1, GE plot 1, SM plot 1…
lfmc
lfmc %>%
  count(leaf.type)

lfmc %>%
  count(time)

lfmc %>%
  count(plot)


lfmc %>%
  count(site)

lfmc %>%
  count(group)

lfmc %>%
  count(site, plot)
lfmc %>%
  ggplot(aes(time, lfmc)) +
  geom_point() +
  facet_wrap(~group)

lfmc.gd <- groupedData(lfmc ~ time | group, data = lfmc, order.groups = FALSE)
lfmc.gd
Grouped Data: lfmc ~ time | group
           leaf.type time plot site  lfmc     group
1        M. spinosum    1    1    P 238.4 SM plot 1
2        M. spinosum    1    1    P 269.8 SM plot 1
3        M. spinosum    1    1    P 325.2 SM plot 1
4   S. bracteolactus    1    1    P 326.7 SS plot 1
5   S. bracteolactus    1    1    P 257.3 SS plot 1
6   S. bracteolactus    1    1    P 279.3 SS plot 1
7            Grass E    1    1    P  91.9 GE plot 1
8            Grass E    1    1    P  95.0 GE plot 1
32       M. spinosum   13    1    P 215.8 SM plot 1
33       M. spinosum   13    1    P 216.7 SM plot 1
34       M. spinosum   13    1    P 205.3 SM plot 1
35  S. bracteolactus   13    1    P 208.5 SS plot 1
36  S. bracteolactus   13    1    P 225.7 SS plot 1
37  S. bracteolactus   13    1    P 246.8 SS plot 1
38           Grass E   13    1    P  81.3 GE plot 1
39           Grass E   13    1    P  82.4 GE plot 1
40           Grass E   13    1    P  78.1 GE plot 1
68       M. spinosum   26    1    P 194.4 SM plot 1
69       M. spinosum   26    1    P 166.0 SM plot 1
70       M. spinosum   26    1    P 212.0 SM plot 1
71  S. bracteolactus   26    1    P 232.2 SS plot 1
72  S. bracteolactus   26    1    P 190.1 SS plot 1
73  S. bracteolactus   26    1    P 241.6 SS plot 1
74           Grass E   26    1    P  53.9 GE plot 1
75           Grass E   26    1    P  71.9 GE plot 1
76           Grass E   26    1    P  63.1 GE plot 1
104      M. spinosum   36    1    P 178.0 SM plot 1
105      M. spinosum   36    1    P 165.9 SM plot 1
106      M. spinosum   36    1    P 144.0 SM plot 1
107 S. bracteolactus   36    1    P 144.1 SS plot 1
108 S. bracteolactus   36    1    P 142.7 SS plot 1
109 S. bracteolactus   36    1    P 170.0 SS plot 1
110          Grass E   36    1    P  45.7 GE plot 1
111          Grass E   36    1    P  37.3 GE plot 1
112          Grass E   36    1    P  33.3 GE plot 1
140      M. spinosum   57    1    P 103.5 SM plot 1
141      M. spinosum   57    1    P 115.2 SM plot 1
142      M. spinosum   57    1    P 115.2 SM plot 1
143 S. bracteolactus   57    1    P  76.7 SS plot 1
144 S. bracteolactus   57    1    P 110.5 SS plot 1
145 S. bracteolactus   57    1    P 148.5 SS plot 1
146          Grass E   57    1    P  23.0 GE plot 1
147          Grass E   57    1    P  23.3 GE plot 1
148          Grass E   57    1    P  25.6 GE plot 1
176      M. spinosum   68    1    P  65.8 SM plot 1
177      M. spinosum   68    1    P  86.9 SM plot 1
178      M. spinosum   68    1    P  79.5 SM plot 1
179 S. bracteolactus   68    1    P  70.3 SS plot 1
180 S. bracteolactus   68    1    P  67.5 SS plot 1
181 S. bracteolactus   68    1    P  87.6 SS plot 1
182          Grass E   68    1    P  14.0 GE plot 1
183          Grass E   68    1    P  24.3 GE plot 1
184          Grass E   68    1    P  17.1 GE plot 1
212      M. spinosum   80    1    P  67.0 SM plot 1
213      M. spinosum   80    1    P  70.5 SM plot 1
214      M. spinosum   80    1    P  57.9 SM plot 1
215 S. bracteolactus   80    1    P  70.6 SS plot 1
216 S. bracteolactus   80    1    P  55.4 SS plot 1
217 S. bracteolactus   80    1    P  86.0 SS plot 1
218          Grass E   80    1    P  11.4 GE plot 1
219          Grass E   80    1    P  16.2 GE plot 1
220          Grass E   80    1    P  12.8 GE plot 1
9        M. spinosum    1    2    P 251.7 SM plot 2
10       M. spinosum    1    2    P 238.7 SM plot 2
11       M. spinosum    1    2    P 252.6 SM plot 2
12  S. bracteolactus    1    2    P 257.6 SS plot 2
13  S. bracteolactus    1    2    P 265.7 SS plot 2
14  S. bracteolactus    1    2    P 273.9 SS plot 2
15           Grass E    1    2    P  60.6 GE plot 2
16           Grass E    1    2    P  53.7 GE plot 2
41       M. spinosum   13    2    P 228.0 SM plot 2
42       M. spinosum   13    2    P 231.4 SM plot 2
43       M. spinosum   13    2    P 226.6 SM plot 2
44  S. bracteolactus   13    2    P 222.7 SS plot 2
45  S. bracteolactus   13    2    P 243.9 SS plot 2
46  S. bracteolactus   13    2    P 218.6 SS plot 2
47           Grass E   13    2    P  50.5 GE plot 2
48           Grass E   13    2    P  64.3 GE plot 2
49           Grass E   13    2    P  64.3 GE plot 2
77       M. spinosum   26    2    P 185.8 SM plot 2
78       M. spinosum   26    2    P 169.7 SM plot 2
79       M. spinosum   26    2    P 154.8 SM plot 2
80  S. bracteolactus   26    2    P 207.2 SS plot 2
81  S. bracteolactus   26    2    P 200.0 SS plot 2
82  S. bracteolactus   26    2    P 204.9 SS plot 2
83           Grass E   26    2    P  50.5 GE plot 2
84           Grass E   26    2    P  43.2 GE plot 2
85           Grass E   26    2    P  52.7 GE plot 2
113      M. spinosum   36    2    P 115.9 SM plot 2
114      M. spinosum   36    2    P 124.4 SM plot 2
115      M. spinosum   36    2    P 130.7 SM plot 2
116 S. bracteolactus   36    2    P 214.9 SS plot 2
117 S. bracteolactus   36    2    P 159.2 SS plot 2
118 S. bracteolactus   36    2    P 135.9 SS plot 2
119          Grass E   36    2    P  32.1 GE plot 2
120          Grass E   36    2    P  34.0 GE plot 2
121          Grass E   36    2    P  39.2 GE plot 2
149      M. spinosum   57    2    P  96.0 SM plot 2
150      M. spinosum   57    2    P  76.7 SM plot 2
151      M. spinosum   57    2    P 101.2 SM plot 2
152 S. bracteolactus   57    2    P  80.6 SS plot 2
153 S. bracteolactus   57    2    P  83.7 SS plot 2
154 S. bracteolactus   57    2    P 108.0 SS plot 2
155          Grass E   57    2    P  20.4 GE plot 2
156          Grass E   57    2    P  23.7 GE plot 2
157          Grass E   57    2    P  21.7 GE plot 2
185      M. spinosum   68    2    P  62.2 SM plot 2
186      M. spinosum   68    2    P  70.6 SM plot 2
187      M. spinosum   68    2    P  76.6 SM plot 2
188 S. bracteolactus   68    2    P  72.8 SS plot 2
189 S. bracteolactus   68    2    P  85.2 SS plot 2
190 S. bracteolactus   68    2    P  75.4 SS plot 2
191          Grass E   68    2    P  15.4 GE plot 2
192          Grass E   68    2    P  11.8 GE plot 2
193          Grass E   68    2    P  15.3 GE plot 2
221      M. spinosum   80    2    P  38.1 SM plot 2
222      M. spinosum   80    2    P  61.2 SM plot 2
223      M. spinosum   80    2    P  53.2 SM plot 2
224 S. bracteolactus   80    2    P  72.2 SS plot 2
225 S. bracteolactus   80    2    P  64.2 SS plot 2
226 S. bracteolactus   80    2    P  61.4 SS plot 2
227          Grass E   80    2    P  13.9 GE plot 2
228          Grass E   80    2    P  11.2 GE plot 2
229          Grass E   80    2    P  13.3 GE plot 2
17       M. spinosum    1    3    P 224.3 SM plot 3
18       M. spinosum    1    3    P 207.3 SM plot 3
19  S. bracteolactus    1    3    P 240.3 SS plot 3
20  S. bracteolactus    1    3    P 215.8 SS plot 3
21           Grass E    1    3    P  75.3 GE plot 3
22           Grass E    1    3    P  65.7 GE plot 3
50       M. spinosum   13    3    P 240.8 SM plot 3
51       M. spinosum   13    3    P 178.4 SM plot 3
52       M. spinosum   13    3    P 231.2 SM plot 3
53  S. bracteolactus   13    3    P 218.6 SS plot 3
54  S. bracteolactus   13    3    P 244.1 SS plot 3
55  S. bracteolactus   13    3    P 204.4 SS plot 3
56           Grass E   13    3    P  69.3 GE plot 3
57           Grass E   13    3    P  77.2 GE plot 3
58           Grass E   13    3    P  62.5 GE plot 3
86       M. spinosum   26    3    P 154.7 SM plot 3
87       M. spinosum   26    3    P 221.0 SM plot 3
88       M. spinosum   26    3    P 231.7 SM plot 3
89  S. bracteolactus   26    3    P 171.4 SS plot 3
90  S. bracteolactus   26    3    P 209.4 SS plot 3
91  S. bracteolactus   26    3    P 180.5 SS plot 3
92           Grass E   26    3    P  54.7 GE plot 3
93           Grass E   26    3    P  45.9 GE plot 3
94           Grass E   26    3    P  41.1 GE plot 3
122      M. spinosum   36    3    P 164.7 SM plot 3
123      M. spinosum   36    3    P 142.9 SM plot 3
124      M. spinosum   36    3    P 161.4 SM plot 3
125 S. bracteolactus   36    3    P 152.0 SS plot 3
126 S. bracteolactus   36    3    P 114.6 SS plot 3
127 S. bracteolactus   36    3    P 183.8 SS plot 3
128          Grass E   36    3    P  24.4 GE plot 3
129          Grass E   36    3    P  28.0 GE plot 3
130          Grass E   36    3    P  26.4 GE plot 3
158      M. spinosum   57    3    P 129.1 SM plot 3
159      M. spinosum   57    3    P 119.9 SM plot 3
160      M. spinosum   57    3    P 148.1 SM plot 3
161 S. bracteolactus   57    3    P 110.6 SS plot 3
162 S. bracteolactus   57    3    P  54.1 SS plot 3
163 S. bracteolactus   57    3    P  95.1 SS plot 3
164          Grass E   57    3    P  31.3 GE plot 3
165          Grass E   57    3    P  22.5 GE plot 3
166          Grass E   57    3    P  27.4 GE plot 3
 [ reached 'max' / getOption("max.print") -- omitted 81 rows ]
lfmc.gd$plot <- with(lfmc.gd, factor(plot, levels = c("4", "5", "6", "1", "2", "3")))
lfmc.gd
Grouped Data: lfmc ~ time | group
           leaf.type time plot site  lfmc     group
1        M. spinosum    1    1    P 238.4 SM plot 1
2        M. spinosum    1    1    P 269.8 SM plot 1
3        M. spinosum    1    1    P 325.2 SM plot 1
4   S. bracteolactus    1    1    P 326.7 SS plot 1
5   S. bracteolactus    1    1    P 257.3 SS plot 1
6   S. bracteolactus    1    1    P 279.3 SS plot 1
7            Grass E    1    1    P  91.9 GE plot 1
8            Grass E    1    1    P  95.0 GE plot 1
32       M. spinosum   13    1    P 215.8 SM plot 1
33       M. spinosum   13    1    P 216.7 SM plot 1
34       M. spinosum   13    1    P 205.3 SM plot 1
35  S. bracteolactus   13    1    P 208.5 SS plot 1
36  S. bracteolactus   13    1    P 225.7 SS plot 1
37  S. bracteolactus   13    1    P 246.8 SS plot 1
38           Grass E   13    1    P  81.3 GE plot 1
39           Grass E   13    1    P  82.4 GE plot 1
40           Grass E   13    1    P  78.1 GE plot 1
68       M. spinosum   26    1    P 194.4 SM plot 1
69       M. spinosum   26    1    P 166.0 SM plot 1
70       M. spinosum   26    1    P 212.0 SM plot 1
71  S. bracteolactus   26    1    P 232.2 SS plot 1
72  S. bracteolactus   26    1    P 190.1 SS plot 1
73  S. bracteolactus   26    1    P 241.6 SS plot 1
74           Grass E   26    1    P  53.9 GE plot 1
75           Grass E   26    1    P  71.9 GE plot 1
76           Grass E   26    1    P  63.1 GE plot 1
104      M. spinosum   36    1    P 178.0 SM plot 1
105      M. spinosum   36    1    P 165.9 SM plot 1
106      M. spinosum   36    1    P 144.0 SM plot 1
107 S. bracteolactus   36    1    P 144.1 SS plot 1
108 S. bracteolactus   36    1    P 142.7 SS plot 1
109 S. bracteolactus   36    1    P 170.0 SS plot 1
110          Grass E   36    1    P  45.7 GE plot 1
111          Grass E   36    1    P  37.3 GE plot 1
112          Grass E   36    1    P  33.3 GE plot 1
140      M. spinosum   57    1    P 103.5 SM plot 1
141      M. spinosum   57    1    P 115.2 SM plot 1
142      M. spinosum   57    1    P 115.2 SM plot 1
143 S. bracteolactus   57    1    P  76.7 SS plot 1
144 S. bracteolactus   57    1    P 110.5 SS plot 1
145 S. bracteolactus   57    1    P 148.5 SS plot 1
146          Grass E   57    1    P  23.0 GE plot 1
147          Grass E   57    1    P  23.3 GE plot 1
148          Grass E   57    1    P  25.6 GE plot 1
176      M. spinosum   68    1    P  65.8 SM plot 1
177      M. spinosum   68    1    P  86.9 SM plot 1
178      M. spinosum   68    1    P  79.5 SM plot 1
179 S. bracteolactus   68    1    P  70.3 SS plot 1
180 S. bracteolactus   68    1    P  67.5 SS plot 1
181 S. bracteolactus   68    1    P  87.6 SS plot 1
182          Grass E   68    1    P  14.0 GE plot 1
183          Grass E   68    1    P  24.3 GE plot 1
184          Grass E   68    1    P  17.1 GE plot 1
212      M. spinosum   80    1    P  67.0 SM plot 1
213      M. spinosum   80    1    P  70.5 SM plot 1
214      M. spinosum   80    1    P  57.9 SM plot 1
215 S. bracteolactus   80    1    P  70.6 SS plot 1
216 S. bracteolactus   80    1    P  55.4 SS plot 1
217 S. bracteolactus   80    1    P  86.0 SS plot 1
218          Grass E   80    1    P  11.4 GE plot 1
219          Grass E   80    1    P  16.2 GE plot 1
220          Grass E   80    1    P  12.8 GE plot 1
9        M. spinosum    1    2    P 251.7 SM plot 2
10       M. spinosum    1    2    P 238.7 SM plot 2
11       M. spinosum    1    2    P 252.6 SM plot 2
12  S. bracteolactus    1    2    P 257.6 SS plot 2
13  S. bracteolactus    1    2    P 265.7 SS plot 2
14  S. bracteolactus    1    2    P 273.9 SS plot 2
15           Grass E    1    2    P  60.6 GE plot 2
16           Grass E    1    2    P  53.7 GE plot 2
41       M. spinosum   13    2    P 228.0 SM plot 2
42       M. spinosum   13    2    P 231.4 SM plot 2
43       M. spinosum   13    2    P 226.6 SM plot 2
44  S. bracteolactus   13    2    P 222.7 SS plot 2
45  S. bracteolactus   13    2    P 243.9 SS plot 2
46  S. bracteolactus   13    2    P 218.6 SS plot 2
47           Grass E   13    2    P  50.5 GE plot 2
48           Grass E   13    2    P  64.3 GE plot 2
49           Grass E   13    2    P  64.3 GE plot 2
77       M. spinosum   26    2    P 185.8 SM plot 2
78       M. spinosum   26    2    P 169.7 SM plot 2
79       M. spinosum   26    2    P 154.8 SM plot 2
80  S. bracteolactus   26    2    P 207.2 SS plot 2
81  S. bracteolactus   26    2    P 200.0 SS plot 2
82  S. bracteolactus   26    2    P 204.9 SS plot 2
83           Grass E   26    2    P  50.5 GE plot 2
84           Grass E   26    2    P  43.2 GE plot 2
85           Grass E   26    2    P  52.7 GE plot 2
113      M. spinosum   36    2    P 115.9 SM plot 2
114      M. spinosum   36    2    P 124.4 SM plot 2
115      M. spinosum   36    2    P 130.7 SM plot 2
116 S. bracteolactus   36    2    P 214.9 SS plot 2
117 S. bracteolactus   36    2    P 159.2 SS plot 2
118 S. bracteolactus   36    2    P 135.9 SS plot 2
119          Grass E   36    2    P  32.1 GE plot 2
120          Grass E   36    2    P  34.0 GE plot 2
121          Grass E   36    2    P  39.2 GE plot 2
149      M. spinosum   57    2    P  96.0 SM plot 2
150      M. spinosum   57    2    P  76.7 SM plot 2
151      M. spinosum   57    2    P 101.2 SM plot 2
152 S. bracteolactus   57    2    P  80.6 SS plot 2
153 S. bracteolactus   57    2    P  83.7 SS plot 2
154 S. bracteolactus   57    2    P 108.0 SS plot 2
155          Grass E   57    2    P  20.4 GE plot 2
156          Grass E   57    2    P  23.7 GE plot 2
157          Grass E   57    2    P  21.7 GE plot 2
185      M. spinosum   68    2    P  62.2 SM plot 2
186      M. spinosum   68    2    P  70.6 SM plot 2
187      M. spinosum   68    2    P  76.6 SM plot 2
188 S. bracteolactus   68    2    P  72.8 SS plot 2
189 S. bracteolactus   68    2    P  85.2 SS plot 2
190 S. bracteolactus   68    2    P  75.4 SS plot 2
191          Grass E   68    2    P  15.4 GE plot 2
192          Grass E   68    2    P  11.8 GE plot 2
193          Grass E   68    2    P  15.3 GE plot 2
221      M. spinosum   80    2    P  38.1 SM plot 2
222      M. spinosum   80    2    P  61.2 SM plot 2
223      M. spinosum   80    2    P  53.2 SM plot 2
224 S. bracteolactus   80    2    P  72.2 SS plot 2
225 S. bracteolactus   80    2    P  64.2 SS plot 2
226 S. bracteolactus   80    2    P  61.4 SS plot 2
227          Grass E   80    2    P  13.9 GE plot 2
228          Grass E   80    2    P  11.2 GE plot 2
229          Grass E   80    2    P  13.3 GE plot 2
17       M. spinosum    1    3    P 224.3 SM plot 3
18       M. spinosum    1    3    P 207.3 SM plot 3
19  S. bracteolactus    1    3    P 240.3 SS plot 3
20  S. bracteolactus    1    3    P 215.8 SS plot 3
21           Grass E    1    3    P  75.3 GE plot 3
22           Grass E    1    3    P  65.7 GE plot 3
50       M. spinosum   13    3    P 240.8 SM plot 3
51       M. spinosum   13    3    P 178.4 SM plot 3
52       M. spinosum   13    3    P 231.2 SM plot 3
53  S. bracteolactus   13    3    P 218.6 SS plot 3
54  S. bracteolactus   13    3    P 244.1 SS plot 3
55  S. bracteolactus   13    3    P 204.4 SS plot 3
56           Grass E   13    3    P  69.3 GE plot 3
57           Grass E   13    3    P  77.2 GE plot 3
58           Grass E   13    3    P  62.5 GE plot 3
86       M. spinosum   26    3    P 154.7 SM plot 3
87       M. spinosum   26    3    P 221.0 SM plot 3
88       M. spinosum   26    3    P 231.7 SM plot 3
89  S. bracteolactus   26    3    P 171.4 SS plot 3
90  S. bracteolactus   26    3    P 209.4 SS plot 3
91  S. bracteolactus   26    3    P 180.5 SS plot 3
92           Grass E   26    3    P  54.7 GE plot 3
93           Grass E   26    3    P  45.9 GE plot 3
94           Grass E   26    3    P  41.1 GE plot 3
122      M. spinosum   36    3    P 164.7 SM plot 3
123      M. spinosum   36    3    P 142.9 SM plot 3
124      M. spinosum   36    3    P 161.4 SM plot 3
125 S. bracteolactus   36    3    P 152.0 SS plot 3
126 S. bracteolactus   36    3    P 114.6 SS plot 3
127 S. bracteolactus   36    3    P 183.8 SS plot 3
128          Grass E   36    3    P  24.4 GE plot 3
129          Grass E   36    3    P  28.0 GE plot 3
130          Grass E   36    3    P  26.4 GE plot 3
158      M. spinosum   57    3    P 129.1 SM plot 3
159      M. spinosum   57    3    P 119.9 SM plot 3
160      M. spinosum   57    3    P 148.1 SM plot 3
161 S. bracteolactus   57    3    P 110.6 SS plot 3
162 S. bracteolactus   57    3    P  54.1 SS plot 3
163 S. bracteolactus   57    3    P  95.1 SS plot 3
164          Grass E   57    3    P  31.3 GE plot 3
165          Grass E   57    3    P  22.5 GE plot 3
166          Grass E   57    3    P  27.4 GE plot 3
 [ reached 'max' / getOption("max.print") -- omitted 81 rows ]
ggplot(lfmc, aes(time, lfmc)) +
  geom_point() +
  geom_smooth(method = "loess") +
  facet_wrap(~leaf.type) +
  xlab("Time (days)") +
  ylab("LMFC (%)")
`geom_smooth()` using formula 'y ~ x'

NONLINEAR MIXED-EFFECTS MODEL (M1)

nlsL <- nlsList(lfmc ~ SSdlf(time, A, w, m, s), data = lfmc.gd)
Warning: 3 errors caught in nls(model, data = data, control = controlvals).  The error messages and their frequencies are

step factor 0.000488281 reduced below 'minFactor' of 0.000976562 
                                                               1 
                     number of iterations exceeded maximum of 50 
                                                               2 
coef(nlsL)
plot(intervals(nlsL))

nl.re <- nlme(nlsL, control = nlmeControl(maxIter = 5000, msMaxIter = 1500))
Warning in (function (model, data = sys.frame(sys.parent()), fixed, random,  :
  Iteration 1, LME step: nlminb() did not converge (code = 1). PORT message: function evaluation limit reached without convergence (9)
nl.re.1 <- update(nl.re, random = pdDiag(A + w + m + s ~ 1))
fxf <- fixef(nl.re.1)
fxf
        A         w         m         s 
193.45267  21.45740  31.48270 -22.59711 
# nl.me <- update(nl.re.1, fixed = list(A + w + m + s ~ leaf.type),
#                 start = c(A = fxf[1],0,0,0,
#                           w = fxf[2],0,0,0,
#                           m = fxf[3],0,0,0,
#                           s = fxf[4],0,0,0))
nl.re.2 <- update(nl.re.1, random = pdDiag(A + w + m ~ 1))
nl.re.2
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -1070.648
  Fixed: list(A ~ 1, w ~ 1, m ~ 1, s ~ 1) 
        A         w         m         s 
193.46656  21.43199  31.48851 -22.60602 

Random effects:
 Formula: list(A ~ 1, w ~ 1, m ~ 1)
 Level: group
 Structure: Diagonal
               A        w        m Residual
StdDev: 119.3792 10.94166 8.927276 15.26574

Number of Observations: 247
Number of Groups: 12 
fxf <- fixef(nl.re.2)
fxf
        A         w         m         s 
193.46656  21.43199  31.48851 -22.60602 
nl.me.1 <- update(
  nl.re.2,
  fixed = list(A + w + m + s ~ leaf.type),
  start = c(A = fxf[1], 0, 0, 0, w = fxf[2], 0, 0, 0, m = fxf[3], 0, 0, 0, s = fxf[4], 0, 0, 0)
)
nl.me.1
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -1034.375
  Fixed: list(A + w + m + s ~ leaf.type) 
              A.(Intercept)          A.leaf.typeGrass E      A.leaf.typeM. spinosum A.leaf.typeS. bracteolactus 
                  23.135414                   53.294908                  296.620596                  263.313118 
              w.(Intercept)          w.leaf.typeGrass E      w.leaf.typeM. spinosum w.leaf.typeS. bracteolactus 
                  50.291577                  -34.313924                  -26.732287                    2.135665 
              m.(Intercept)          m.leaf.typeGrass E      m.leaf.typeM. spinosum m.leaf.typeS. bracteolactus 
                  51.533983                  -21.920404                  -20.647403                  -18.170359 
              s.(Intercept)          s.leaf.typeGrass E      s.leaf.typeM. spinosum s.leaf.typeS. bracteolactus 
                  17.605794                  -25.926031                  -43.529450                  -33.666226 

Random effects:
 Formula: list(A ~ 1, w ~ 1, m ~ 1)
 Level: group
 Structure: Diagonal
        A.(Intercept) w.(Intercept) m.(Intercept) Residual
StdDev:      8.112813   0.001236897      2.468334 15.38607

Number of Observations: 247
Number of Groups: 12 
AIC(nl.re.2, nl.me.1)
IC_tab(nl.re.2, nl.me.1)

The small StdDev of “w” suggests to remove the random-effect of plot on this parameter.

nl.me.2 <- update(nl.me.1, random = pdDiag(A + m ~ 1))
IC_tab(nl.me.1, nl.me.2)

in addition, if the random-effects on “m” is also removed.

nl.me.3 <- update(nl.me.2, random = pdDiag(A ~ 1))
nl.me.3
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -1035.19
  Fixed: list(A + w + m + s ~ leaf.type) 
              A.(Intercept)          A.leaf.typeGrass E      A.leaf.typeM. spinosum A.leaf.typeS. bracteolactus 
                  23.071500                   53.950821                  408.256907                  258.698322 
              w.(Intercept)          w.leaf.typeGrass E      w.leaf.typeM. spinosum w.leaf.typeS. bracteolactus 
                  50.311106                  -34.625804                  -79.299174                    3.717378 
              m.(Intercept)          m.leaf.typeGrass E      m.leaf.typeM. spinosum m.leaf.typeS. bracteolactus 
                  51.608221                  -21.975396                  -32.109785                  -17.665025 
              s.(Intercept)          s.leaf.typeGrass E      s.leaf.typeM. spinosum s.leaf.typeS. bracteolactus 
                  17.683236                  -26.518320                  -60.147105                  -33.113726 

Random effects:
 Formula: A ~ 1 | group
        A.(Intercept) Residual
StdDev:       10.6035 15.53457

Number of Observations: 247
Number of Groups: 12 

the model fit is similar, supporting the simpler model.

IC_tab(nl.me.2, nl.me.3)
plot(nl.me.3)

hist(residuals(nl.me.3, type = "normalized"), breaks = 20)

rse <- residuals(nl.me.3) / sigma(nl.me.3)
par(pty = "s")
qqnorm(rse)
qqline(rse)

Heterocedasticity: variance modeling —-

nl.me.3.vm <- update(nl.me.3, weights = varIdent(form = ~ 1 | leaf.type))
nl.me.3.vm
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -965.5628
  Fixed: list(A + w + m + s ~ leaf.type) 
              A.(Intercept)          A.leaf.typeGrass E      A.leaf.typeM. spinosum A.leaf.typeS. bracteolactus 
                 50.1522851                  28.9553861                 346.5816334                 231.1569667 
              w.(Intercept)          w.leaf.typeGrass E      w.leaf.typeM. spinosum w.leaf.typeS. bracteolactus 
                 24.6693607                  -9.4837876                 -35.8327942                  29.7533014 
              m.(Intercept)          m.leaf.typeGrass E      m.leaf.typeM. spinosum m.leaf.typeS. bracteolactus 
                 48.9613646                 -19.9780200                 -26.2920794                 -14.9909275 
              s.(Intercept)          s.leaf.typeGrass E      s.leaf.typeM. spinosum s.leaf.typeS. bracteolactus 
                -16.2294127                   6.3520529                 -21.0664455                   0.9389323 

Random effects:
 Formula: A ~ 1 | group
        A.(Intercept) Residual
StdDev:      9.412263 6.824811

Variance function:
 Structure: Different standard deviations per stratum
 Formula: ~1 | leaf.type 
 Parameter estimates:
         Grass W          Grass E      M. spinosum S. bracteolactus 
       1.0000000        0.8547835        3.1456077        3.0262116 
Number of Observations: 247
Number of Groups: 12 
IC_tab(nl.me.3, nl.me.3.vm)

여기서 가중잔차는 weight를 이용해서 계산해야 제대로 plot이 됨

res <-
  tibble(
    fitted = fitted(nl.me.3.vm),
    resid = resid(nl.me.3.vm), # 일반 잔차
    nresid = resid(nl.me.3.vm, type = "normalized"), # 표준 잔차
    weight = diag(var_cov(nl.me.3.vm)), # 가중치
    wresid = sqrt(1 / weight) * resid # 가중잔차 수동계산
  ) %>%
  bind_cols(lfmc.gd)
res
res %>%
  ggplot(aes(fitted, nresid)) +
  geom_point(aes(color = leaf.type))

res %>%
  ggplot(aes(fitted, wresid)) +
  geom_point(aes(color = leaf.type))

plot(nl.me.3.vm)

hist(residuals(nl.me.3.vm, type = "normalized"), main = "", xlab = "Residuals")

the residuals look fine now.

qqnorm(resid(nl.me.3.vm, type = "normalized"))
qqline(resid(nl.me.3.vm, type = "normalized"))

Temporal correlation check. This plot suggests there is no need for modeling the correlation of the residuals.

plot(ACF(nl.me.3.vm), alpha = 0.01)

Evaluation of Fixed Effects

the model is fitted by Maximum likelihood

nl.me.ml <- update(nl.me.3.vm, method = "ML")
nl.me.ml
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -965.5628
  Fixed: list(A + w + m + s ~ leaf.type) 
              A.(Intercept)          A.leaf.typeGrass E      A.leaf.typeM. spinosum A.leaf.typeS. bracteolactus 
                 50.1522851                  28.9553861                 346.5816334                 231.1569667 
              w.(Intercept)          w.leaf.typeGrass E      w.leaf.typeM. spinosum w.leaf.typeS. bracteolactus 
                 24.6693607                  -9.4837876                 -35.8327942                  29.7533014 
              m.(Intercept)          m.leaf.typeGrass E      m.leaf.typeM. spinosum m.leaf.typeS. bracteolactus 
                 48.9613646                 -19.9780200                 -26.2920794                 -14.9909275 
              s.(Intercept)          s.leaf.typeGrass E      s.leaf.typeM. spinosum s.leaf.typeS. bracteolactus 
                -16.2294127                   6.3520529                 -21.0664455                   0.9389323 

Random effects:
 Formula: A ~ 1 | group
        A.(Intercept) Residual
StdDev:      9.412263 6.824811

Variance function:
 Structure: Different standard deviations per stratum
 Formula: ~1 | leaf.type 
 Parameter estimates:
         Grass W          Grass E      M. spinosum S. bracteolactus 
       1.0000000        0.8547835        3.1456077        3.0262116 
Number of Observations: 247
Number of Groups: 12 

the model is reduced removing the fixed-effect of “leaf type” on “s”

nl.me.ml.1 <- update(nl.me.ml,
  fixed = list(A + w + m ~ leaf.type, s ~ 1),
  start = c(A = fxf[1], 0, 0, 0, w = fxf[2], 0, 0, 0, m = fxf[3], 0, 0, 0, s = fxf[4])
)
nl.me.ml.1
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -969.2617
  Fixed: list(A + w + m ~ leaf.type, s ~ 1) 
              A.(Intercept)          A.leaf.typeGrass E      A.leaf.typeM. spinosum A.leaf.typeS. bracteolactus 
                   49.94562                    39.60979                   219.45652                   232.11909 
              w.(Intercept)          w.leaf.typeGrass E      w.leaf.typeM. spinosum w.leaf.typeS. bracteolactus 
                   25.24501                   -14.34293                    35.36122                    28.70389 
              m.(Intercept)          m.leaf.typeGrass E      m.leaf.typeM. spinosum m.leaf.typeS. bracteolactus 
                   48.30390                   -21.59477                   -15.57765                   -14.37454 
                          s 
                  -15.44291 

Random effects:
 Formula: A ~ 1 | group
        A.(Intercept) Residual
StdDev:      9.612472 6.825604

Variance function:
 Structure: Different standard deviations per stratum
 Formula: ~1 | leaf.type 
 Parameter estimates:
         Grass W          Grass E      M. spinosum S. bracteolactus 
       1.0000000        0.8714808        3.2750696        3.0238816 
Number of Observations: 247
Number of Groups: 12 

the AIC comparison suggests that the more complex model is not justified (fits are similar)

IC_tab(nl.me.ml, nl.me.ml.1)

the model is reduced removing the fixed-effect of “leaf type” on “m”

nl.me.ml.2 <- update(nl.me.ml,
  fixed = list(A + w ~ leaf.type, m + s ~ 1),
  start = c(A = fxf[1], 0, 0, 0, w = fxf[2], 0, 0, 0, m = fxf[3], s = fxf[4])
)
nl.me.ml.2
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -972.7465
  Fixed: list(A + w ~ leaf.type, m + s ~ 1) 
              A.(Intercept)          A.leaf.typeGrass E      A.leaf.typeM. spinosum A.leaf.typeS. bracteolactus 
                   54.25730                    30.93396                   223.27914                   240.31542 
              w.(Intercept)          w.leaf.typeGrass E      w.leaf.typeM. spinosum w.leaf.typeS. bracteolactus 
                   29.05803                   -20.68663                    31.68646                    26.99003 
                          m                           s 
                   30.91655                   -16.15435 

Random effects:
 Formula: A ~ 1 | group
        A.(Intercept) Residual
StdDev:      9.216113 7.028345

Variance function:
 Structure: Different standard deviations per stratum
 Formula: ~1 | leaf.type 
 Parameter estimates:
         Grass W          Grass E      M. spinosum S. bracteolactus 
       1.0000000        0.8820414        3.1693973        2.9643303 
Number of Observations: 247
Number of Groups: 12 

the model is reduced removing the fixed-effect of “leaf type” on “W”.

nl.me.ml.3 <- update(nl.me.ml,
  fixed = list(A ~ leaf.type, w + m + s ~ 1),
  start = c(A = fxf[1], 0, 0, 0, w = fxf[2], m = fxf[3], s = fxf[4])
)
nl.me.ml.3
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -1017.91
  Fixed: list(A ~ leaf.type, w + m + s ~ 1) 
              A.(Intercept)          A.leaf.typeGrass E      A.leaf.typeM. spinosum A.leaf.typeS. bracteolactus 
                   63.92200                    14.27252                   251.99435                   265.42764 
                          w                           m                           s 
                   17.33019                    32.65673                   -25.16298 

Random effects:
 Formula: A ~ 1 | group
        A.(Intercept) Residual
StdDev:      8.138947 7.795662

Variance function:
 Structure: Different standard deviations per stratum
 Formula: ~1 | leaf.type 
 Parameter estimates:
         Grass W          Grass E      M. spinosum S. bracteolactus 
        1.000000         1.584887         2.813985         2.731415 
Number of Observations: 247
Number of Groups: 12 
nl.me.ml.4 <- update(nl.me.ml,
  fixed = list(w ~ leaf.type, A + m + s ~ 1),
  start = c(A = fxf[1], w = fxf[2], 0, 0, 0, m = fxf[3], s = fxf[4])
)
nl.me.ml.4
Nonlinear mixed-effects model fit by maximum likelihood
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
  Log-likelihood: -999.8224
  Fixed: list(w ~ leaf.type, A + m + s ~ 1) 
              w.(Intercept)          w.leaf.typeGrass E      w.leaf.typeM. spinosum w.leaf.typeS. bracteolactus 
                   29.16029                   -20.38028                    33.19430                    28.72481 
                          A                           m                           s 
                  175.78131                    31.08598                   -15.70297 

Random effects:
 Formula: A ~ 1 | group
               A Residual
StdDev: 108.0655 7.042308

Variance function:
 Structure: Different standard deviations per stratum
 Formula: ~1 | leaf.type 
 Parameter estimates:
         Grass W          Grass E      M. spinosum S. bracteolactus 
       1.0000000        0.8770721        3.1889180        2.9304615 
Number of Observations: 247
Number of Groups: 12 
IC_tab(nl.me.ml.1, nl.me.ml.2, nl.me.ml.3, nl.me.ml.4)

Final Model

This model assumes:

  • “A” and “w” vary with leaf type (fixed-effects).
  • “A” varies with plot (random-effects).
  • “m” and “s” are unique for all the plots and leaf types.
  • Residual variance depends on “leaf type”
  • There is no temporal dependence
  • Residuals following a normal distribution
M1 <- update(nl.me.ml.2, method = "REML")
summary(M1)
Nonlinear mixed-effects model fit by REML
  Model: lfmc ~ SSdlf(time, A, w, m, s) 
  Data: lfmc.gd 
       AIC      BIC   logLik
  1931.668 1983.689 -950.834

Random effects:
 Formula: A ~ 1 | group
        A.(Intercept) Residual
StdDev:      9.408521 7.162899

Variance function:
 Structure: Different standard deviations per stratum
 Formula: ~1 | leaf.type 
 Parameter estimates:
         Grass W          Grass E      M. spinosum S. bracteolactus 
       1.0000000        0.8851764        3.1796102        2.9647393 
Fixed effects:  list(A + w ~ leaf.type, m + s ~ 1) 
                                Value Std.Error  DF   t-value p-value
A.(Intercept)                54.25350  5.988181 226  9.060097   0e+00
A.leaf.typeGrass E           30.92293  8.835210 226  3.499966   6e-04
A.leaf.typeM. spinosum      223.24504 15.915557 226 14.026844   0e+00
A.leaf.typeS. bracteolactus 240.27654 16.857356 226 14.253513   0e+00
w.(Intercept)                29.05581  1.712960 226 16.962334   0e+00
w.leaf.typeGrass E          -20.68938  2.592469 226 -7.980569   0e+00
w.leaf.typeM. spinosum       31.67297  7.753176 226  4.085161   1e-04
w.leaf.typeS. bracteolactus  26.97508  8.071111 226  3.342177   1e-03
m                            30.92881  2.065196 226 14.976213   0e+00
s                           -16.15406  2.365392 226 -6.829339   0e+00
 Correlation: 
                            A.(In) A.l.GE A..M.s A..S.b w.(In) w.l.GE w..M.s w..S.b m     
A.leaf.typeGrass E          -0.527                                                        
A.leaf.typeM. spinosum      -0.146  0.535                                                 
A.leaf.typeS. bracteolactus -0.115  0.537  0.746                                          
w.(Intercept)               -0.217 -0.034 -0.200 -0.216                                   
w.leaf.typeGrass E          -0.035 -0.283 -0.382 -0.399 -0.258                            
w.leaf.typeM. spinosum      -0.118 -0.227 -0.547 -0.454  0.157  0.573                     
w.leaf.typeS. bracteolactus -0.129 -0.241 -0.462 -0.575  0.188  0.601  0.640              
m                           -0.217 -0.324 -0.633 -0.666  0.092  0.145  0.161  0.172       
s                           -0.246 -0.354 -0.710 -0.748  0.416  0.575  0.702  0.752  0.544

Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-2.27999661 -0.59657156 -0.02042043  0.57273326  3.12018001 

Number of Observations: 247
Number of Groups: 12 
fixef(M1)
              A.(Intercept)          A.leaf.typeGrass E      A.leaf.typeM. spinosum A.leaf.typeS. bracteolactus 
                   54.25350                    30.92293                   223.24504                   240.27654 
              w.(Intercept)          w.leaf.typeGrass E      w.leaf.typeM. spinosum w.leaf.typeS. bracteolactus 
                   29.05581                   -20.68938                    31.67297                    26.97508 
                          m                           s 
                   30.92881                   -16.15406 
ranef(M1)
          A.(Intercept)
GW plot 4     -1.915443
GW plot 5     -1.923198
GW plot 6      3.838641
GE plot 1     15.652459
SM plot 1      6.967303
SS plot 1      9.166925
GE plot 2    -11.066197
SM plot 2     -5.358910
SS plot 2      2.442781
GE plot 3     -4.586263
SM plot 3     -1.608393
SS plot 3    -11.609707
intervals(M1)
Approximate 95% confidence intervals

 Fixed effects:
                                lower      est.     upper
A.(Intercept)                42.45370  54.25350  66.05331
A.leaf.typeGrass E           13.51301  30.92293  48.33286
A.leaf.typeM. spinosum      191.88318 223.24504 254.60691
A.leaf.typeS. bracteolactus 207.05884 240.27654 273.49423
w.(Intercept)                25.68039  29.05581  32.43122
w.leaf.typeGrass E          -25.79788 -20.68938 -15.58088
w.leaf.typeM. spinosum       16.39521  31.67297  46.95073
w.leaf.typeS. bracteolactus  11.07083  26.97508  42.87934
m                            26.85931  30.92881  34.99831
s                           -20.81510 -16.15406 -11.49302

 Random Effects:
  Level: group 
                   lower     est.    upper
sd(A.(Intercept)) 5.0352 9.408521 17.58029

 Variance function:
                    lower      est.    upper
Grass E          0.675418 0.8851764 1.160077
M. spinosum      2.451780 3.1796102 4.123502
S. bracteolactus 2.285695 2.9647393 3.845517

 Within-group standard error:
   lower     est.    upper 
5.966409 7.162899 8.599330 

Part II

# 2.2.1 - Overall predictions (Table 3) ----
Agw <- fixef(M1)[1] # "A" for grasses in the W site (GW)
Age <- fixef(M1)[1] + fixef(M1)[2] # "A" for grasses in the E site (GE)
Asm <- fixef(M1)[1] + fixef(M1)[3] # "A" for Mullinum spinosum (SM)
Ass <- fixef(M1)[1] + fixef(M1)[4] # "A" for Senecio filaginoides (SS)
Wgw <- fixef(M1)[5] # "w" for grasses in the W site (GW)
Wge <- fixef(M1)[5] + fixef(M1)[6] # "w" for grasses in the E site (GE)
Wsm <- fixef(M1)[5] + fixef(M1)[7] # "w" for Mullinum spinosum (SM)
Wss <- fixef(M1)[5] + fixef(M1)[8] # "w" for Senecio filaginoides (SS)
M <- fixef(M1)[9] # "m" for all the leaf types
S <- fixef(M1)[10] # "s" for all the leaf types

GW <- subset(lfmc, leaf.type == "Grass W")
GE <- subset(lfmc, leaf.type == "Grass E")
SM <- subset(lfmc, leaf.type == "M. spinosum")
SS <- subset(lfmc, leaf.type == "S. bracteolactus")

par(mfrow = c(2, 2), cex = 0.75)
# Grasses W
plot(GW$time, GW$lfmc, xlab = "", ylab = "LFMC (%)", ylim = c(0, 100), main = "Grasses W site", col = "darkgreen")
curve(((Agw - Wgw) / (1 + exp((M - x) / S)) + Wgw), lwd = 2, lty = 1, add = T, col = "darkgreen")
# Grasses E
plot(GE$time, GE$lfmc, xlab = "", ylab = "LFMC (%)", ylim = c(0, 100), main = "Grasses E site", col = "green")
curve(((Age - Wge) / (1 + exp((M - x) / S)) + Wge), lwd = 2, lty = 1, add = T, col = "green")
# M. Spinosum
plot(SM$time, SM$lfmc, xlab = "Time (days)", ylab = "LFMC (%)", ylim = c(0, 350), main = "M. spinosum", font.main = 4, col = "darkorange")
curve(((Asm - Wsm) / (1 + exp((M - x) / S)) + Wsm), lwd = 2, lty = 1, add = T, col = "darkorange")
# S. filaginoides
plot(SS$time, SS$lfmc, xlab = "Time (days)", ylab = "LFMC (%)", ylim = c(0, 350), main = "S. bracteolactus", font.main = 4, col = "darkred")
curve(((Ass - Wss) / (1 + exp((M - x) / S)) + Wss), lwd = 2, lty = 1, add = T, col = "darkred")

ranef(M1)
          A.(Intercept)
GW plot 4     -1.915443
GW plot 5     -1.923198
GW plot 6      3.838641
GE plot 1     15.652459
SM plot 1      6.967303
SS plot 1      9.166925
GE plot 2    -11.066197
SM plot 2     -5.358910
SS plot 2      2.442781
GE plot 3     -4.586263
SM plot 3     -1.608393
SS plot 3    -11.609707
# 2.2.2 - Predictions at plot level (Table 3) ----
Agw.p4 <- fixef(M1)[1] + ranef(M1)[1, 1] # "A" for grasses in the plot 4 (GW)
Agw.p5 <- fixef(M1)[1] + ranef(M1)[2, 1] # "A" for grasses in the plot 5 (GW)
Agw.p6 <- fixef(M1)[1] + ranef(M1)[3, 1] # "A" for grasses in the plot 6 (GW)
Age.p1 <- fixef(M1)[1] + fixef(M1)[2] + ranef(M1)[4, 1] # "A" for grasses in the plot 1 (GW)
Age.p2 <- fixef(M1)[1] + fixef(M1)[2] + ranef(M1)[7, 1] # "A" for grasses in the plot 2 (GW)
Age.p3 <- fixef(M1)[1] + fixef(M1)[2] + ranef(M1)[10, 1] # "A" for grasses in the plot 3 (GW)
Asm.p1 <- fixef(M1)[1] + fixef(M1)[3] + ranef(M1)[5, 1] # "A" for M. spinosum in the plot 1 (SM)
Asm.p2 <- fixef(M1)[1] + fixef(M1)[3] + ranef(M1)[8, 1] # "A" for M. spinosum in the plot 2 (SM)
Asm.p3 <- fixef(M1)[1] + fixef(M1)[3] + ranef(M1)[11, 1] # "A" for M. spinosum in the plot 3 (SM)
Ass.p1 <- fixef(M1)[1] + fixef(M1)[4] + ranef(M1)[6, 1] # "A" for S. filaginoides in the plot 1 (SM)
Ass.p2 <- fixef(M1)[1] + fixef(M1)[4] + ranef(M1)[9, 1] # "A" for S. filaginoides in the plot 2 (SM)
Ass.p3 <- fixef(M1)[1] + fixef(M1)[4] + ranef(M1)[12, 1] # "A" for S. filaginoides in the plot 3 (SM)

par(mfrow = c(2, 2), cex = 0.75) # (Figure 4)
# Grasses W site
plot(GW$time, GW$lfmc, xlab = "", ylab = "LFMC (%)", ylim = c(0, 100), main = "Grasses W site", col = "darkgreen")
curve(((Agw.p4 - Wgw) / (1 + exp((M - x) / S)) + Wgw), lwd = 1, lty = 2, add = T, col = "darkgreen") # plot 4
curve(((Agw.p5 - Wgw) / (1 + exp((M - x) / S)) + Wgw), lwd = 1, lty = 2, add = T, col = "darkgreen") # plot 5
curve(((Agw.p6 - Wgw) / (1 + exp((M - x) / S)) + Wgw), lwd = 1, lty = 2, add = T, col = "darkgreen") # plot 6
# Grasses E site
plot(GE$time, GE$lfmc, xlab = "", ylab = "LFMC (%)", ylim = c(0, 100), main = "Grasses E site", col = "green")
curve(((Age.p1 - Wge) / (1 + exp((M - x) / S)) + Wge), lwd = 1, lty = 2, add = T, col = "green") # plot 1
curve(((Age.p2 - Wge) / (1 + exp((M - x) / S)) + Wge), lwd = 1, lty = 2, add = T, col = "green") # plot 2
curve(((Age.p3 - Wge) / (1 + exp((M - x) / S)) + Wge), lwd = 1, lty = 2, add = T, col = "green") # plot 3
# M. Spinosum (E site)
plot(SM$time, SM$lfmc, xlab = "Time (days)", ylab = "LFMC (%)", ylim = c(0, 350), main = "M. spinosum", font.main = 4, col = "darkorange")
curve(((Asm.p1 - Wsm) / (1 + exp((M - x) / S)) + Wsm), lwd = 1, lty = 2, add = T, col = "darkorange") # plot 1
curve(((Asm.p2 - Wsm) / (1 + exp((M - x) / S)) + Wsm), lwd = 1, lty = 2, add = T, col = "darkorange") # plot 2
curve(((Asm.p3 - Wsm) / (1 + exp((M - x) / S)) + Wsm), lwd = 1, lty = 2, add = T, col = "darkorange") # plot 3
# S. filaginoides (E site)
plot(SS$time, SS$lfmc, xlab = "Time (days)", ylab = "LFMC (%)", ylim = c(0, 350), main = "S. bracteolactus", font.main = 4, col = "darkred")
curve(((Ass.p1 - Wss) / (1 + exp((M - x) / S)) + Wss), lwd = 1, lty = 2, add = T, col = "darkred") # plot 1
curve(((Ass.p2 - Wss) / (1 + exp((M - x) / S)) + Wss), lwd = 1, lty = 2, add = T, col = "darkred") # plot 2
curve(((Ass.p3 - Wss) / (1 + exp((M - x) / S)) + Wss), lwd = 1, lty = 2, add = T, col = "darkred") # plot 3

# 2.2.3 - Derivatives (drying speed) ----
lfmc.GW <- expression((Agw - Wgw) / (1 + exp((M - x) / S)) + Wgw) # LFMC dynamics for grasses in the W site
ds.GW <- deriv(lfmc.GW, "x", function.arg = TRUE) # Drying speed for grasses in the W site

lfmc.GE <- expression((Age - Wge) / (1 + exp((M - x) / S)) + Wge) # LFMC dynamics for grasses in the E site
ds.GE <- deriv(lfmc.GE, "x", function.arg = TRUE) # Drying speed for grasses in the E site

lfmc.SM <- expression((Asm - Wsm) / (1 + exp((M - x) / S)) + Wsm) # LFMC dynamics for M. spinosum
ds.SM <- deriv(lfmc.SM, "x", function.arg = TRUE) # Drying speed for M. spinosum

lfmc.SS <- expression((Ass - Wss) / (1 + exp((M - x) / S)) + Wss) # LFMC dynamics for S. filaginoides
ds.SS <- deriv(lfmc.SS, "x", function.arg = TRUE) # # Drying speed for S. filaginoides

xvec <- seq(0, 100, length = 1000)
y.ds.GW <- ds.GW(xvec)
y.ds.GE <- ds.GE(xvec)
y.ds.SM <- ds.SM(xvec)
y.ds.SS <- ds.SS(xvec)
par(mfrow = c(2, 2), cex = 0.75) # (Figure 4)
plot(xvec, y.ds.GW,
  xlim = c(0, 80), ylim = c(0, 4),
  ylab = "Drying speed", xlab = "Time (days)", type = "n", main = "Grasses in the W site"
)
lines(xvec, -1 * (attr(y.ds.GW, "grad")), lwd = 2, lty = 1, col = "darkgreen")
plot(xvec, y.ds.GE,
  xlim = c(0, 80), ylim = c(0, 4),
  ylab = "Drying speed", xlab = "Time (days)", type = "n", main = "Grasses in the E site"
)
lines(xvec, -1 * (attr(y.ds.GE, "grad")), lwd = 2, lty = 1, col = "green")
plot(xvec, y.ds.SM,
  xlim = c(0, 80), ylim = c(0, 4),
  ylab = "Drying speed", xlab = "Time (days)", type = "n", main = "M. spinosum", font.main = 4
)
lines(xvec, -1 * (attr(y.ds.SM, "grad")), lwd = 2, lty = 1, col = "darkorange")
plot(xvec, y.ds.SS,
  xlim = c(0, 80), ylim = c(0, 4),
  ylab = "Drying speed", xlab = "Time (days)", type = "n", main = "S. filaginoides", font.main = 4
)
lines(xvec, -1 * (attr(y.ds.SS, "grad")), lwd = 2, lty = 1, col = "darkred")

# 2.3 ---- Testing differences among leaf types -------------------------------------------

# 2.3.1 - Parameter "A" ----
contrast(emmeans(M1, ~leaf.type, param = "A"), "pairwise")
 contrast                       estimate    SE  df t.ratio p.value
 Grass W - Grass E                 -30.9  8.84 226  -3.500  0.0031
 Grass W - M. spinosum            -223.2 15.92 226 -14.027  <.0001
 Grass W - S. bracteolactus       -240.3 16.86 226 -14.254  <.0001
 Grass E - M. spinosum            -192.3 13.45 226 -14.294  <.0001
 Grass E - S. bracteolactus       -209.4 14.22 226 -14.718  <.0001
 M. spinosum - S. bracteolactus    -17.0 11.71 226  -1.454  0.4671

P value adjustment: tukey method for comparing a family of 4 estimates 
# 2.3.2 - Parameter "w" ----
contrast(emmeans(M1, ~leaf.type, param = "w"), "pairwise")
 contrast                       estimate   SE  df t.ratio p.value
 Grass W - Grass E                  20.7 2.59 226   7.981  <.0001
 Grass W - M. spinosum             -31.7 7.75 226  -4.085  0.0004
 Grass W - S. bracteolactus        -27.0 8.07 226  -3.342  0.0053
 Grass E - M. spinosum             -52.4 6.62 226  -7.914  <.0001
 Grass E - S. bracteolactus        -47.7 6.84 226  -6.974  <.0001
 M. spinosum - S. bracteolactus      4.7 6.72 226   0.699  0.8974

P value adjustment: tukey method for comparing a family of 4 estimates 

3) NONLINEAR FIXED-EFFEC

Response function: lfmc = (A - w) / (1 + exp((m - time)/s))) + w

# 3.1.1 - Base nonlinear model ----
nl.fe.0 <- nls(lfmc ~ ((A - w) / (1 + exp((m - time) / s))) + w,
  start = c(A = 15, m = 30, s = -17, w = 5),
  data = lfmc
)
# here, time is the only predictor variable
nl.fe.0
Nonlinear regression model
  model: lfmc ~ ((A - w)/(1 + exp((m - time)/s))) + w
   data: lfmc
     A      m      s      w 
194.97  29.66 -21.17  27.93 
 residual sum-of-squares: 1029759

Number of iterations to convergence: 8 
Achieved convergence tolerance: 9.57e-07
b <- coef(nl.fe.0) # coefficients of the base model
b
        A         m         s         w 
194.97252  29.66004 -21.17131  27.92798 
# 3.1.2 - Leaf type fixed effect ----
nl.fe.1 <- nls(lfmc ~ ((A[leaf.type] - w[leaf.type]) / (1 + exp((m - time) / s))) + w[leaf.type],
  start = list(A = rep(b[1], 4), m = 25, s = -10, w = rep(b[4], 4)),
  data = lfmc
) # and fit is made using the base model's coefficients as start values
nl.fe.1
Nonlinear regression model
  model: lfmc ~ ((A[leaf.type] - w[leaf.type])/(1 + exp((m - time)/s))) +     w[leaf.type]
   data: lfmc
     A1      A2      A3      A4       m       s      w1      w2      w3      w4 
 56.726  92.483 298.835 318.207  30.431 -20.912  27.172   2.795  44.528  38.298 
 residual sum-of-squares: 67786

Number of iterations to convergence: 5 
Achieved convergence tolerance: 5.958e-06
b1 <- coef(nl.fe.1) # coefficients of the model with leaf type and time as predictor variables
b1
        A1         A2         A3         A4          m          s         w1         w2         w3         w4 
 56.725716  92.482865 298.834954 318.207067  30.430747 -20.911576  27.171667   2.795385  44.527958  38.297817 
# 3.1.3 - Plot fixed effect ----
nl.fe.2 <- nls(lfmc ~ ((A[group] - w[group]) / (1 + exp((m - time) / s))) + w[group],
  start = list(A = rep(c(b1[1], b1[2], b1[3], b1[4]), 3), m = 25, s = -10, w = rep(c(b1[7], b1[8], b1[9], b1[10]), 3)),
  data = lfmc
) # the model is fit using "b1" as the start values
nl.fe.2
Nonlinear regression model
  model: lfmc ~ ((A[group] - w[group])/(1 + exp((m - time)/s))) + w[group]
   data: lfmc
     A1      A2      A3      A4      A5      A6      A7      A8      A9     A10     A11     A12       m       s 
 53.683  54.282  61.906 111.562 314.131 333.831  77.465 298.161 321.433  87.553 279.437 292.443  30.413 -20.696 
     w1      w2      w3      w4      w5      w6      w7      w8      w9     w10     w11     w12 
 28.369  27.480  25.962   0.890  43.524  41.223   6.091  26.610  39.495   2.265  66.661  38.222 
 residual sum-of-squares: 53914

Number of iterations to convergence: 6 
Achieved convergence tolerance: 7.249e-06
IC_tab(nl.fe.0, nl.fe.1, nl.fe.2)
# Residual checking:
par(mfrow = c(1, 1))
plot(nl.fe.2) # heterocedasticity in the residuals

hist(resid(nl.fe.2), main = "", xlab = "Residuals") # slightly skew distribution

qqnorm(resid(nl.fe.2))
qqline(resid(nl.fe.2))

This model assumes:

  • “A” and “w” vary with leaf type and plot (fixed-effects).
  • “m” and “s” are unique for all the plots and leaf types.
  • Variance is homogeneous
  • There is no temporal dependence
  • Residuals following a normal distribution
M2 <- nl.fe.2
summary(M2)

Formula: lfmc ~ ((A[group] - w[group])/(1 + exp((m - time)/s))) + w[group]

Parameters:
    Estimate Std. Error t value Pr(>|t|)    
A1    53.683      8.602   6.240 2.20e-09 ***
A2    54.282      8.636   6.285 1.72e-09 ***
A3    61.906      8.881   6.971 3.59e-11 ***
A4   111.561     13.291   8.394 5.65e-15 ***
A5   314.131     24.840  12.646  < 2e-16 ***
A6   333.831     26.640  12.531  < 2e-16 ***
A7    77.465     10.997   7.044 2.34e-11 ***
A8   298.161     24.917  11.966  < 2e-16 ***
A9   321.433     25.765  12.475  < 2e-16 ***
A10   87.553     11.744   7.455 2.01e-12 ***
A11  279.437     20.829  13.416  < 2e-16 ***
A12  292.443     24.178  12.095  < 2e-16 ***
m     30.413      2.896  10.504  < 2e-16 ***
s    -20.695      3.606  -5.739 3.11e-08 ***
w1    28.369      6.537   4.340 2.17e-05 ***
w2    27.480      6.548   4.197 3.93e-05 ***
w3    25.962      6.633   3.914 0.000121 ***
w4     0.890      8.173   0.109 0.913389    
w5    43.524     13.567   3.208 0.001535 ** 
w6    41.223     14.429   2.857 0.004686 ** 
w7     6.091      7.262   0.839 0.402467    
w8    26.610     13.604   1.956 0.051725 .  
w9    39.495     14.009   2.819 0.005252 ** 
w10    2.265      7.551   0.300 0.764455    
w11   66.661     11.478   5.808 2.19e-08 ***
w12   38.222     13.031   2.933 0.003708 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 15.62 on 221 degrees of freedom

Number of iterations to convergence: 6 
Achieved convergence tolerance: 7.249e-06

4 - LINEAR MIXED-EFFECTS MODEL (M3)

Response function: lfmc = β0 + β1xTime + β2xGE + β3xSM + β4xSS β5xGExTime + β6xSMxTime + β7xSSxTime

where GE, SM, and SS are dummy variables (0 or 1) created to indicate differences between the base level of “leaf type”, i.e., grasses from the W site [GW], and the remaining levels.

4.1.1 Base linear mixed-effect mode

# 4.1.1 - Base linear mixed-effect model ----
l.me <- lme(lfmc ~ time * leaf.type, random = ~ 1 | plot, data = lfmc)

# Residual checking:
plot(l.me) # heterocedasticity in the residuals

plot(ACF(l.me, resType = "normalized", na.action = na.omit), alpha = 0.05, grid = TRUE) # temporal correlation is observed at lag 1

# 4.1.2 - Variance modelling ----
l.me.vm <- update(l.me, weights = varIdent(form = ~ 1 | leaf.type))

# Residual checking:
plot(l.me.vm) # variances there seem to be well modeled

AIC(l.me, l.me.vm) # modeling the variance improves the model fit
# 4.1.3 - Temporal correlation modelling ----
l.me.vm.tm <- update(l.me.vm, correlation = corARMA(p = 1, q = 1, form = ~1))
plot(ACF(l.me.vm.tm, resType = "normalized", na.action = na.omit), alpha = 0.05, grid = TRUE) # the temporal dependence is modeled

AIC(l.me.vm, l.me.vm.tm) # modeling the temporal correlation improves the model fit
plot(l.me.vm.tm)

plot(ACF(l.me.vm.tm, resType = "normalized", na.action = na.omit), alpha = 0.05, grid = TRUE)

hist(resid(l.me.vm.tm, type = "normalized"), main = "", xlab = "Residuals")

qqnorm(resid(l.me.vm.tm, type = "normalized"))
qqline(resid(l.me.vm.tm, type = "normalized"))

# 4.1.4 - Evaluation of the fixed-effects ----
l.me.ml <- update(l.me.vm.tm, method = "ML") # the model is fitted by Maximum likelihood
l.me.ml.1 <- update(l.me.ml, ~ . - time:leaf.type) # the model is reduced removing the interaction "time x leaf type"
IC_tab(l.me.ml, l.me.ml.1) # the AIC comparison suggests the interaction term to be important
M3 <- l.me.vm.tm
summary(M3)
Linear mixed-effects model fit by REML
  Data: lfmc 
       AIC     BIC    logLik
  1998.483 2050.63 -984.2416

Random effects:
 Formula: ~1 | plot
        (Intercept) Residual
StdDev:    3.750488 7.842661

Correlation Structure: ARMA(1,1)
 Formula: ~1 | plot 
 Parameter estimate(s):
      Phi1     Theta1 
 0.9076001 -0.7567060 
Variance function:
 Structure: Different standard deviations per stratum
 Formula: ~1 | leaf.type 
 Parameter estimates:
         Grass W      M. spinosum S. bracteolactus          Grass E 
        1.000000         2.963571         3.123740         1.070897 
Fixed effects:  lfmc ~ time * leaf.type 
                                   Value Std.Error  DF    t-value p-value
(Intercept)                     50.14877  3.478211 234  14.417977       0
time                            -0.27535  0.048114 234  -5.722819       0
leaf.typeGrass E                24.40614  4.908714 234   4.972002       0
leaf.typeM. spinosum           194.85292  8.539072 234  22.818980       0
leaf.typeS. bracteolactus      209.78398  8.838810 234  23.734415       0
time:leaf.typeGrass E           -0.56367  0.071704 234  -7.861042       0
time:leaf.typeM. spinosum       -2.02767  0.152372 234 -13.307359       0
time:leaf.typeS. bracteolactus  -2.29338  0.160937 234 -14.250193       0
 Correlation: 
                               (Intr) time   lf.tGE lf.M.s lf.S.b tm:.GE t:.M.s
time                           -0.557                                          
leaf.typeGrass E               -0.709  0.395                                   
leaf.typeM. spinosum           -0.407  0.227  0.654                            
leaf.typeS. bracteolactus      -0.394  0.219  0.653  0.666                     
time:leaf.typeGrass E           0.374 -0.671 -0.592 -0.430 -0.429              
time:leaf.typeM. spinosum       0.176 -0.316 -0.299 -0.742 -0.409  0.546       
time:leaf.typeS. bracteolactus  0.166 -0.299 -0.318 -0.441 -0.743  0.562  0.567

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max 
-2.3449313 -0.6438200 -0.1011963  0.6158869  3.3469739 

Number of Observations: 247
Number of Groups: 6 
M4 <- lm(lfmc ~ time * group, data = lfmc)
summary(M4)

Call:
lm(formula = lfmc ~ time * group, data = lfmc)

Residuals:
    Min      1Q  Median      3Q     Max 
-49.673  -7.042   0.466   7.002  66.332 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)          48.79275    6.44309   7.573 9.60e-13 ***
time                 -0.24478    0.13305  -1.840  0.06714 .  
groupGW plot 5        0.49817    9.11190   0.055  0.95645    
groupGW plot 6        6.14984    9.11190   0.675  0.50042    
groupGE plot 1       40.32675    9.49655   4.246 3.19e-05 ***
groupSM plot 1      212.68278    9.11190  23.341  < 2e-16 ***
groupSS plot 1      227.42468    9.11190  24.959  < 2e-16 ***
groupGE plot 2       14.39850    9.49655   1.516  0.13089    
groupSM plot 2      195.54780    9.11190  21.461  < 2e-16 ***
groupSS plot 2      217.12857    9.11190  23.829  < 2e-16 ***
groupGE plot 3       21.68513    9.49655   2.283  0.02334 *  
groupSM plot 3      189.79613    9.49655  19.986  < 2e-16 ***
groupSS plot 3      192.53490    9.49655  20.274  < 2e-16 ***
time:groupGW plot 5  -0.01905    0.18816  -0.101  0.91944    
time:groupGW plot 6  -0.10231    0.18816  -0.544  0.58718    
time:groupGE plot 1  -0.80129    0.19357  -4.139 4.94e-05 ***
time:groupSM plot 1  -2.36256    0.18816 -12.556  < 2e-16 ***
time:groupSS plot 1  -2.55767    0.18816 -13.593  < 2e-16 ***
time:groupGE plot 2  -0.43459    0.19357  -2.245  0.02574 *  
time:groupSM plot 2  -2.34722    0.18816 -12.474  < 2e-16 ***
time:groupSS plot 2  -2.45552    0.18816 -13.050  < 2e-16 ***
time:groupGE plot 3  -0.56657    0.19357  -2.927  0.00378 ** 
time:groupSM plot 3  -1.82099    0.19357  -9.407  < 2e-16 ***
time:groupSS plot 3  -2.16847    0.19357 -11.202  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 16.51 on 223 degrees of freedom
Multiple R-squared:  0.9586,    Adjusted R-squared:  0.9544 
F-statistic: 224.7 on 23 and 223 DF,  p-value: < 2.2e-16
# Residual checking:
par(mfrow = c(2, 2))
plot(M4) # a clear pattern in the residuals is observed (model assumptions are not met)

# 6 - NULL MODEL (M5) ######################################################################################

M5 <- lm(lfmc ~ 1, data = lfmc)
summary(M5)

Call:
lm(formula = lfmc ~ 1, data = lfmc)

Residuals:
   Min     1Q Median     3Q    Max 
-88.75 -58.05 -31.45  52.66 231.06 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   95.645      4.919   19.45   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 77.3 on 246 degrees of freedom
IC_tab(M2, M3, M4, M5) # , weights = T, delta = TRUE, base=T, sort = TRUE) # the nonlinear mixed-effects model is clearly the best
IC_tab(M1)
par(mfrow = c(2, 2), cex = 0.65) # (Figure 5)
v1 <- c(1, 2, 3, 4)
v2 <- c("GW", "GE", "SM", "SS")

# Nonlinear mixed-effects model (M1)
rM1 <- resid(M1, type = "normalized")
fM1 <- fitted(M1)
plot(fM1, rM1, xlab = "Fitted LFMC (%)", ylab = "Residuals", ylim = c(-3, 3))
abline(a = 0, b = 0, lwd = 2)
boxplot(rM1 ~ lfmc$leaf.type, ylab = "Residuals", xlab = "Leaf type", xaxt = "n", ylim = c(-3, 3))
axis(side = 1, at = v1, labels = v2, las = 1)
plot(rM1 ~ lfmc$time, ylab = "Residuals", xlab = "Time (days)", ylim = c(-3, 3))
abline(a = 0, b = 0, lw = 2)
qqnorm(rM1, main = "", xlab = "Theoretical quantiles", ylab = "Sample quantiles")
qqline(rM1)

# Nonlinear fixed-effects model (M2)
par(mfrow = c(2, 2), cex = 0.65) # (Figure 5)

rM2 <- resid(M2)
fM2 <- fitted(M2)
plot(fM2, rM2, xlab = "Fitted LFMC (%)", ylab = "Residuals")
abline(a = 0, b = 0, lwd = 2)
boxplot(rM2 ~ lfmc$leaf.type, ylab = "Residuals", xlab = "Leaf type", xaxt = "n")
axis(side = 1, at = v1, labels = v2, las = 1)
plot(rM2 ~ lfmc$time, ylab = "Residuals", xlab = "Time (days)")
abline(a = 0, b = 0, lw = 2)
qqnorm(rM2, main = "", xlab = "Theoretical quantiles", ylab = "Sample quantiles")
qqline(rM2)

par(mfrow = c(2, 2), cex = 0.65)

# Linear mixed-effects model (M3)
rM3 <- resid(M3, type = "normalized")
fM3 <- fitted(M3)
plot(fM3, rM3, xlab = "Fitted LFMC (%)", ylab = "Residuals", ylim = c(-3, 3))
abline(a = 0, b = 0, lwd = 2)
boxplot(rM3 ~ lfmc$leaf.type, ylab = "Residuals", xlab = "Leaf type", xaxt = "n", ylim = c(-3, 3))
axis(side = 1, at = v1, labels = v2, las = 1)
plot(rM3 ~ lfmc$time, ylab = "Residuals", xlab = "Time (days)", ylim = c(-3, 3))
abline(a = 0, b = 0, lw = 2)
qqnorm(rM3, main = "", ylab = "Sample quantiles", xlab = "Theoretical quantiles")
qqline(rM3)

par(mfrow = c(2, 2), cex = 0.65)

# Clasical regression (M4)
rM4 <- resid(M4)
fM4 <- fitted(M4)
plot(fM4, rM4, xlab = "Fitted LFMC (%)", ylab = "Residuals")
abline(a = 0, b = 0, lwd = 2)
boxplot(rM4 ~ lfmc$leaf.type, ylab = "Residuals", xlab = "", xaxt = "n")
axis(side = 1, at = v1, labels = v2, las = 1)
plot(rM4 ~ lfmc$time, ylab = "Residuals", xlab = "Time (days)")
abline(a = 0, b = 0, lw = 2)
qqnorm(rM4, main = "", xlab = "Theoretical quantiles", ylab = "Sample quantiles")
qqline(rM4)

—-

LS0tCnRpdGxlOiAiTm9ubGluZWFyIFJlZ3Jlc3Npb24gKEFyY2hvbnRvdWxpcyBhbmQgTWlndWV6KSBwYXBlciIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKPGh0dHBzOi8vZmVtaWd1ZXouZ2l0aHViLmlvL25scmFhLWRvY3MvbmxyYWEtT2RkaS1MRk1DLmh0bWw+Cgo8aHR0cHM6Ly9vbmxpbmVsaWJyYXJ5LndpbGV5LmNvbS9kb2kvcGRmZGlyZWN0LzEwLjEwMDIvZWNlMy41NTQzPgoKYGBge3J9CiMgaW5zdGFsbC5wYWNrYWdlcygiZW1tZWFucyIpCgpsaWJyYXJ5KG5sbWUpCmxpYnJhcnkobmxyYWEpCmxpYnJhcnkoZ2dwbG90MikKbGlicmFyeShlbW1lYW5zKSAjIGVtbWVhbnM6OmNvbnRyYXN0KCksIGVtbWVhbnM6OmVtbWVhbnMoKQpsaWJyYXJ5KHRpZHl2ZXJzZSkKYGBgCgojIFN0ZXBzIGluIGZpdHRpbmcgYSBub25saW5lYXIgbWl4ZWQgbW9kZWwKCiMjIFBhcnQgSS4KCmBgYHtyfQpkYXRhKGxmbWMpCmdsaW1wc2UobGZtYykKYGBgCgpgYGB7cn0KbGZtYwpgYGAKCmBgYHtyfQpsZm1jICU+JQogIGNvdW50KGxlYWYudHlwZSkKCmxmbWMgJT4lCiAgY291bnQodGltZSkKCmxmbWMgJT4lCiAgY291bnQocGxvdCkKCgpsZm1jICU+JQogIGNvdW50KHNpdGUpCgpsZm1jICU+JQogIGNvdW50KGdyb3VwKQoKbGZtYyAlPiUKICBjb3VudChzaXRlLCBwbG90KQpgYGAKCmBgYHtyfQpsZm1jICU+JQogIGdncGxvdChhZXModGltZSwgbGZtYykpICsKICBnZW9tX3BvaW50KCkgKwogIGZhY2V0X3dyYXAofmdyb3VwKQpgYGAKCmBgYHtyfQpsZm1jLmdkIDwtIGdyb3VwZWREYXRhKGxmbWMgfiB0aW1lIHwgZ3JvdXAsIGRhdGEgPSBsZm1jLCBvcmRlci5ncm91cHMgPSBGQUxTRSkKbGZtYy5nZApgYGAKCmBgYHtyfQpsZm1jLmdkJHBsb3QgPC0gd2l0aChsZm1jLmdkLCBmYWN0b3IocGxvdCwgbGV2ZWxzID0gYygiNCIsICI1IiwgIjYiLCAiMSIsICIyIiwgIjMiKSkpCmxmbWMuZ2QKYGBgCgpgYGB7cn0KZ2dwbG90KGxmbWMsIGFlcyh0aW1lLCBsZm1jKSkgKwogIGdlb21fcG9pbnQoKSArCiAgZ2VvbV9zbW9vdGgobWV0aG9kID0gImxvZXNzIikgKwogIGZhY2V0X3dyYXAofmxlYWYudHlwZSkgKwogIHhsYWIoIlRpbWUgKGRheXMpIikgKwogIHlsYWIoIkxNRkMgKCUpIikKYGBgCgojIyBOT05MSU5FQVIgTUlYRUQtRUZGRUNUUyBNT0RFTCAoTTEpCgpgYGB7cn0KbmxzTCA8LSBubHNMaXN0KGxmbWMgfiBTU2RsZih0aW1lLCBBLCB3LCBtLCBzKSwgZGF0YSA9IGxmbWMuZ2QpCmNvZWYobmxzTCkKYGBgCgpgYGB7ciwgZmlnLndpZHRoID0gNCwgaGVpZ2h0ID0gM30KcGxvdChpbnRlcnZhbHMobmxzTCkpCmBgYAoKYGBge3J9Cm5sLnJlIDwtIG5sbWUobmxzTCwgY29udHJvbCA9IG5sbWVDb250cm9sKG1heEl0ZXIgPSA1MDAwLCBtc01heEl0ZXIgPSAxNTAwKSkKYGBgCgpgYGB7cn0KbmwucmUuMSA8LSB1cGRhdGUobmwucmUsIHJhbmRvbSA9IHBkRGlhZyhBICsgdyArIG0gKyBzIH4gMSkpCmZ4ZiA8LSBmaXhlZihubC5yZS4xKQpmeGYKYGBgCgpgYGB7cn0KIyBubC5tZSA8LSB1cGRhdGUobmwucmUuMSwgZml4ZWQgPSBsaXN0KEEgKyB3ICsgbSArIHMgfiBsZWFmLnR5cGUpLAojICAgICAgICAgICAgICAgICBzdGFydCA9IGMoQSA9IGZ4ZlsxXSwwLDAsMCwKIyAgICAgICAgICAgICAgICAgICAgICAgICAgIHcgPSBmeGZbMl0sMCwwLDAsCiMgICAgICAgICAgICAgICAgICAgICAgICAgICBtID0gZnhmWzNdLDAsMCwwLAojICAgICAgICAgICAgICAgICAgICAgICAgICAgcyA9IGZ4Zls0XSwwLDAsMCkpCmBgYAoKYGBge3J9Cm5sLnJlLjIgPC0gdXBkYXRlKG5sLnJlLjEsIHJhbmRvbSA9IHBkRGlhZyhBICsgdyArIG0gfiAxKSkKbmwucmUuMgpgYGAKCmBgYHtyfQpmeGYgPC0gZml4ZWYobmwucmUuMikKZnhmCmBgYAoKYGBge3J9Cm5sLm1lLjEgPC0gdXBkYXRlKAogIG5sLnJlLjIsCiAgZml4ZWQgPSBsaXN0KEEgKyB3ICsgbSArIHMgfiBsZWFmLnR5cGUpLAogIHN0YXJ0ID0gYyhBID0gZnhmWzFdLCAwLCAwLCAwLCB3ID0gZnhmWzJdLCAwLCAwLCAwLCBtID0gZnhmWzNdLCAwLCAwLCAwLCBzID0gZnhmWzRdLCAwLCAwLCAwKQopCm5sLm1lLjEKYGBgCgpgYGB7cn0KQUlDKG5sLnJlLjIsIG5sLm1lLjEpCmBgYAoKYGBge3J9CklDX3RhYihubC5yZS4yLCBubC5tZS4xKQpgYGAKClRoZSBzbWFsbCBTdGREZXYgb2YgInciIHN1Z2dlc3RzIHRvIHJlbW92ZSB0aGUgcmFuZG9tLWVmZmVjdCBvZiBwbG90IG9uIHRoaXMgcGFyYW1ldGVyLgoKYGBge3J9Cm5sLm1lLjIgPC0gdXBkYXRlKG5sLm1lLjEsIHJhbmRvbSA9IHBkRGlhZyhBICsgbSB+IDEpKQpJQ190YWIobmwubWUuMSwgbmwubWUuMikKYGBgCgppbiBhZGRpdGlvbiwgaWYgdGhlIHJhbmRvbS1lZmZlY3RzIG9uICJtIiBpcyBhbHNvIHJlbW92ZWQuCgpgYGB7cn0KbmwubWUuMyA8LSB1cGRhdGUobmwubWUuMiwgcmFuZG9tID0gcGREaWFnKEEgfiAxKSkKbmwubWUuMwpgYGAKCnRoZSBtb2RlbCBmaXQgaXMgc2ltaWxhciwgc3VwcG9ydGluZyB0aGUgc2ltcGxlciBtb2RlbC4KCmBgYHtyfQpJQ190YWIobmwubWUuMiwgbmwubWUuMykKYGBgCgpgYGB7cn0KcGxvdChubC5tZS4zKQpgYGAKCmBgYHtyfQpoaXN0KHJlc2lkdWFscyhubC5tZS4zLCB0eXBlID0gIm5vcm1hbGl6ZWQiKSwgYnJlYWtzID0gMjApCmBgYAoKYGBge3J9CnJzZSA8LSByZXNpZHVhbHMobmwubWUuMykgLyBzaWdtYShubC5tZS4zKQpwYXIocHR5ID0gInMiKQpxcW5vcm0ocnNlKQpxcWxpbmUocnNlKQpgYGAKCkhldGVyb2NlZGFzdGljaXR5OiB2YXJpYW5jZSBtb2RlbGluZyAtLS0tCgpgYGB7cn0KbmwubWUuMy52bSA8LSB1cGRhdGUobmwubWUuMywgd2VpZ2h0cyA9IHZhcklkZW50KGZvcm0gPSB+IDEgfCBsZWFmLnR5cGUpKQpubC5tZS4zLnZtCmBgYAoKYGBge3J9CklDX3RhYihubC5tZS4zLCBubC5tZS4zLnZtKQpgYGAKCuyXrOq4sOyEnCDqsIDspJHsnpTssKjripQgd2VpZ2h066W8IOydtOyaqe2VtOyEnCDqs4TsgrDtlbTslbwg7KCc64yA66GcIHBsb3TsnbQg65CoCgpgYGB7cn0KcmVzIDwtCiAgdGliYmxlKAogICAgZml0dGVkID0gZml0dGVkKG5sLm1lLjMudm0pLAogICAgcmVzaWQgPSByZXNpZChubC5tZS4zLnZtKSwgIyDsnbzrsJgg7J6U7LCoCiAgICBucmVzaWQgPSByZXNpZChubC5tZS4zLnZtLCB0eXBlID0gIm5vcm1hbGl6ZWQiKSwgIyDtkZzspIAg7J6U7LCoCiAgICB3ZWlnaHQgPSBkaWFnKHZhcl9jb3YobmwubWUuMy52bSkpLCAjIOqwgOykkey5mAogICAgd3Jlc2lkID0gc3FydCgxIC8gd2VpZ2h0KSAqIHJlc2lkICMg6rCA7KSR7J6U7LCoIOyImOuPmeqzhOyCsAogICkgJT4lCiAgYmluZF9jb2xzKGxmbWMuZ2QpCnJlcwpgYGAKCmBgYHtyfQpyZXMgJT4lCiAgZ2dwbG90KGFlcyhmaXR0ZWQsIG5yZXNpZCkpICsKICBnZW9tX3BvaW50KGFlcyhjb2xvciA9IGxlYWYudHlwZSkpCnJlcyAlPiUKICBnZ3Bsb3QoYWVzKGZpdHRlZCwgd3Jlc2lkKSkgKwogIGdlb21fcG9pbnQoYWVzKGNvbG9yID0gbGVhZi50eXBlKSkKcGxvdChubC5tZS4zLnZtKQpgYGAKCmBgYHtyfQpoaXN0KHJlc2lkdWFscyhubC5tZS4zLnZtLCB0eXBlID0gIm5vcm1hbGl6ZWQiKSwgbWFpbiA9ICIiLCB4bGFiID0gIlJlc2lkdWFscyIpCmBgYAoKdGhlIHJlc2lkdWFscyBsb29rIGZpbmUgbm93LgoKYGBge3J9CnFxbm9ybShyZXNpZChubC5tZS4zLnZtLCB0eXBlID0gIm5vcm1hbGl6ZWQiKSkKcXFsaW5lKHJlc2lkKG5sLm1lLjMudm0sIHR5cGUgPSAibm9ybWFsaXplZCIpKQpgYGAKClRlbXBvcmFsIGNvcnJlbGF0aW9uIGNoZWNrLiBUaGlzIHBsb3Qgc3VnZ2VzdHMgdGhlcmUgaXMgbm8gbmVlZCBmb3IgbW9kZWxpbmcgdGhlIGNvcnJlbGF0aW9uIG9mIHRoZSByZXNpZHVhbHMuCgpgYGB7cn0KcGxvdChBQ0YobmwubWUuMy52bSksIGFscGhhID0gMC4wMSkKYGBgCgojIyMgRXZhbHVhdGlvbiBvZiBGaXhlZCBFZmZlY3RzCgojIyB0aGUgbW9kZWwgaXMgZml0dGVkIGJ5IE1heGltdW0gbGlrZWxpaG9vZAoKYGBge3J9Cm5sLm1lLm1sIDwtIHVwZGF0ZShubC5tZS4zLnZtLCBtZXRob2QgPSAiTUwiKQpubC5tZS5tbApgYGAKCnRoZSBtb2RlbCBpcyByZWR1Y2VkIHJlbW92aW5nIHRoZSBmaXhlZC1lZmZlY3Qgb2YgImxlYWYgdHlwZSIgb24gInMiCgpgYGB7cn0KbmwubWUubWwuMSA8LSB1cGRhdGUobmwubWUubWwsCiAgZml4ZWQgPSBsaXN0KEEgKyB3ICsgbSB+IGxlYWYudHlwZSwgcyB+IDEpLAogIHN0YXJ0ID0gYyhBID0gZnhmWzFdLCAwLCAwLCAwLCB3ID0gZnhmWzJdLCAwLCAwLCAwLCBtID0gZnhmWzNdLCAwLCAwLCAwLCBzID0gZnhmWzRdKQopCm5sLm1lLm1sLjEKYGBgCgp0aGUgQUlDIGNvbXBhcmlzb24gc3VnZ2VzdHMgdGhhdCB0aGUgbW9yZSBjb21wbGV4IG1vZGVsIGlzIG5vdCBqdXN0aWZpZWQgKGZpdHMgYXJlIHNpbWlsYXIpCgpgYGB7cn0KSUNfdGFiKG5sLm1lLm1sLCBubC5tZS5tbC4xKQpgYGAKCnRoZSBtb2RlbCBpcyByZWR1Y2VkIHJlbW92aW5nIHRoZSBmaXhlZC1lZmZlY3Qgb2YgImxlYWYgdHlwZSIgb24gIm0iCgpgYGB7cn0KbmwubWUubWwuMiA8LSB1cGRhdGUobmwubWUubWwsCiAgZml4ZWQgPSBsaXN0KEEgKyB3IH4gbGVhZi50eXBlLCBtICsgcyB+IDEpLAogIHN0YXJ0ID0gYyhBID0gZnhmWzFdLCAwLCAwLCAwLCB3ID0gZnhmWzJdLCAwLCAwLCAwLCBtID0gZnhmWzNdLCBzID0gZnhmWzRdKQopCm5sLm1lLm1sLjIKYGBgCgp0aGUgbW9kZWwgaXMgcmVkdWNlZCByZW1vdmluZyB0aGUgZml4ZWQtZWZmZWN0IG9mICJsZWFmIHR5cGUiIG9uICJXIi4KCmBgYHtyfQpubC5tZS5tbC4zIDwtIHVwZGF0ZShubC5tZS5tbCwKICBmaXhlZCA9IGxpc3QoQSB+IGxlYWYudHlwZSwgdyArIG0gKyBzIH4gMSksCiAgc3RhcnQgPSBjKEEgPSBmeGZbMV0sIDAsIDAsIDAsIHcgPSBmeGZbMl0sIG0gPSBmeGZbM10sIHMgPSBmeGZbNF0pCikKbmwubWUubWwuMwpgYGAKCmBgYHtyfQpubC5tZS5tbC40IDwtIHVwZGF0ZShubC5tZS5tbCwKICBmaXhlZCA9IGxpc3QodyB+IGxlYWYudHlwZSwgQSArIG0gKyBzIH4gMSksCiAgc3RhcnQgPSBjKEEgPSBmeGZbMV0sIHcgPSBmeGZbMl0sIDAsIDAsIDAsIG0gPSBmeGZbM10sIHMgPSBmeGZbNF0pCikKbmwubWUubWwuNApgYGAKCmBgYHtyfQpJQ190YWIobmwubWUubWwuMSwgbmwubWUubWwuMiwgbmwubWUubWwuMywgbmwubWUubWwuNCkKYGBgCgojIyMgRmluYWwgTW9kZWwKClRoaXMgbW9kZWwgYXNzdW1lczoKCi0gICAiQSIgYW5kICJ3IiB2YXJ5IHdpdGggbGVhZiB0eXBlIChmaXhlZC1lZmZlY3RzKS4KLSAgICJBIiB2YXJpZXMgd2l0aCBwbG90IChyYW5kb20tZWZmZWN0cykuCi0gICAibSIgYW5kICJzIiBhcmUgdW5pcXVlIGZvciBhbGwgdGhlIHBsb3RzIGFuZCBsZWFmIHR5cGVzLgotICAgUmVzaWR1YWwgdmFyaWFuY2UgZGVwZW5kcyBvbiAibGVhZiB0eXBlIgotICAgVGhlcmUgaXMgbm8gdGVtcG9yYWwgZGVwZW5kZW5jZQotICAgUmVzaWR1YWxzIGZvbGxvd2luZyBhIG5vcm1hbCBkaXN0cmlidXRpb24KCmBgYHtyLCBwYWdlZC5wcmludCA9IEZBTFNFfQpNMSA8LSB1cGRhdGUobmwubWUubWwuMiwgbWV0aG9kID0gIlJFTUwiKQpzdW1tYXJ5KE0xKQpgYGAKCmBgYHtyfQpmaXhlZihNMSkKYGBgCgpgYGB7cn0KcmFuZWYoTTEpCmBgYAoKYGBge3J9CmludGVydmFscyhNMSkKYGBgCgojIFBhcnQgSUkKCmBgYHtyfQojIDIuMi4xIC0gT3ZlcmFsbCBwcmVkaWN0aW9ucyAoVGFibGUgMykgLS0tLQpBZ3cgPC0gZml4ZWYoTTEpWzFdICMgIkEiIGZvciBncmFzc2VzIGluIHRoZSBXIHNpdGUgKEdXKQpBZ2UgPC0gZml4ZWYoTTEpWzFdICsgZml4ZWYoTTEpWzJdICMgIkEiIGZvciBncmFzc2VzIGluIHRoZSBFIHNpdGUgKEdFKQpBc20gPC0gZml4ZWYoTTEpWzFdICsgZml4ZWYoTTEpWzNdICMgIkEiIGZvciBNdWxsaW51bSBzcGlub3N1bSAoU00pCkFzcyA8LSBmaXhlZihNMSlbMV0gKyBmaXhlZihNMSlbNF0gIyAiQSIgZm9yIFNlbmVjaW8gZmlsYWdpbm9pZGVzIChTUykKV2d3IDwtIGZpeGVmKE0xKVs1XSAjICJ3IiBmb3IgZ3Jhc3NlcyBpbiB0aGUgVyBzaXRlIChHVykKV2dlIDwtIGZpeGVmKE0xKVs1XSArIGZpeGVmKE0xKVs2XSAjICJ3IiBmb3IgZ3Jhc3NlcyBpbiB0aGUgRSBzaXRlIChHRSkKV3NtIDwtIGZpeGVmKE0xKVs1XSArIGZpeGVmKE0xKVs3XSAjICJ3IiBmb3IgTXVsbGludW0gc3Bpbm9zdW0gKFNNKQpXc3MgPC0gZml4ZWYoTTEpWzVdICsgZml4ZWYoTTEpWzhdICMgInciIGZvciBTZW5lY2lvIGZpbGFnaW5vaWRlcyAoU1MpCk0gPC0gZml4ZWYoTTEpWzldICMgIm0iIGZvciBhbGwgdGhlIGxlYWYgdHlwZXMKUyA8LSBmaXhlZihNMSlbMTBdICMgInMiIGZvciBhbGwgdGhlIGxlYWYgdHlwZXMKCkdXIDwtIHN1YnNldChsZm1jLCBsZWFmLnR5cGUgPT0gIkdyYXNzIFciKQpHRSA8LSBzdWJzZXQobGZtYywgbGVhZi50eXBlID09ICJHcmFzcyBFIikKU00gPC0gc3Vic2V0KGxmbWMsIGxlYWYudHlwZSA9PSAiTS4gc3Bpbm9zdW0iKQpTUyA8LSBzdWJzZXQobGZtYywgbGVhZi50eXBlID09ICJTLiBicmFjdGVvbGFjdHVzIikKCnBhcihtZnJvdyA9IGMoMiwgMiksIGNleCA9IDAuNzUpCiMgR3Jhc3NlcyBXCnBsb3QoR1ckdGltZSwgR1ckbGZtYywgeGxhYiA9ICIiLCB5bGFiID0gIkxGTUMgKCUpIiwgeWxpbSA9IGMoMCwgMTAwKSwgbWFpbiA9ICJHcmFzc2VzIFcgc2l0ZSIsIGNvbCA9ICJkYXJrZ3JlZW4iKQpjdXJ2ZSgoKEFndyAtIFdndykgLyAoMSArIGV4cCgoTSAtIHgpIC8gUykpICsgV2d3KSwgbHdkID0gMiwgbHR5ID0gMSwgYWRkID0gVCwgY29sID0gImRhcmtncmVlbiIpCiMgR3Jhc3NlcyBFCnBsb3QoR0UkdGltZSwgR0UkbGZtYywgeGxhYiA9ICIiLCB5bGFiID0gIkxGTUMgKCUpIiwgeWxpbSA9IGMoMCwgMTAwKSwgbWFpbiA9ICJHcmFzc2VzIEUgc2l0ZSIsIGNvbCA9ICJncmVlbiIpCmN1cnZlKCgoQWdlIC0gV2dlKSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXZ2UpLCBsd2QgPSAyLCBsdHkgPSAxLCBhZGQgPSBULCBjb2wgPSAiZ3JlZW4iKQojIE0uIFNwaW5vc3VtCnBsb3QoU00kdGltZSwgU00kbGZtYywgeGxhYiA9ICJUaW1lIChkYXlzKSIsIHlsYWIgPSAiTEZNQyAoJSkiLCB5bGltID0gYygwLCAzNTApLCBtYWluID0gIk0uIHNwaW5vc3VtIiwgZm9udC5tYWluID0gNCwgY29sID0gImRhcmtvcmFuZ2UiKQpjdXJ2ZSgoKEFzbSAtIFdzbSkgLyAoMSArIGV4cCgoTSAtIHgpIC8gUykpICsgV3NtKSwgbHdkID0gMiwgbHR5ID0gMSwgYWRkID0gVCwgY29sID0gImRhcmtvcmFuZ2UiKQojIFMuIGZpbGFnaW5vaWRlcwpwbG90KFNTJHRpbWUsIFNTJGxmbWMsIHhsYWIgPSAiVGltZSAoZGF5cykiLCB5bGFiID0gIkxGTUMgKCUpIiwgeWxpbSA9IGMoMCwgMzUwKSwgbWFpbiA9ICJTLiBicmFjdGVvbGFjdHVzIiwgZm9udC5tYWluID0gNCwgY29sID0gImRhcmtyZWQiKQpjdXJ2ZSgoKEFzcyAtIFdzcykgLyAoMSArIGV4cCgoTSAtIHgpIC8gUykpICsgV3NzKSwgbHdkID0gMiwgbHR5ID0gMSwgYWRkID0gVCwgY29sID0gImRhcmtyZWQiKQpgYGAKCmBgYHtyfQpyYW5lZihNMSkKYGBgCgpgYGB7cn0KIyAyLjIuMiAtIFByZWRpY3Rpb25zIGF0IHBsb3QgbGV2ZWwgKFRhYmxlIDMpIC0tLS0KQWd3LnA0IDwtIGZpeGVmKE0xKVsxXSArIHJhbmVmKE0xKVsxLCAxXSAjICJBIiBmb3IgZ3Jhc3NlcyBpbiB0aGUgcGxvdCA0IChHVykKQWd3LnA1IDwtIGZpeGVmKE0xKVsxXSArIHJhbmVmKE0xKVsyLCAxXSAjICJBIiBmb3IgZ3Jhc3NlcyBpbiB0aGUgcGxvdCA1IChHVykKQWd3LnA2IDwtIGZpeGVmKE0xKVsxXSArIHJhbmVmKE0xKVszLCAxXSAjICJBIiBmb3IgZ3Jhc3NlcyBpbiB0aGUgcGxvdCA2IChHVykKQWdlLnAxIDwtIGZpeGVmKE0xKVsxXSArIGZpeGVmKE0xKVsyXSArIHJhbmVmKE0xKVs0LCAxXSAjICJBIiBmb3IgZ3Jhc3NlcyBpbiB0aGUgcGxvdCAxIChHVykKQWdlLnAyIDwtIGZpeGVmKE0xKVsxXSArIGZpeGVmKE0xKVsyXSArIHJhbmVmKE0xKVs3LCAxXSAjICJBIiBmb3IgZ3Jhc3NlcyBpbiB0aGUgcGxvdCAyIChHVykKQWdlLnAzIDwtIGZpeGVmKE0xKVsxXSArIGZpeGVmKE0xKVsyXSArIHJhbmVmKE0xKVsxMCwgMV0gIyAiQSIgZm9yIGdyYXNzZXMgaW4gdGhlIHBsb3QgMyAoR1cpCkFzbS5wMSA8LSBmaXhlZihNMSlbMV0gKyBmaXhlZihNMSlbM10gKyByYW5lZihNMSlbNSwgMV0gIyAiQSIgZm9yIE0uIHNwaW5vc3VtIGluIHRoZSBwbG90IDEgKFNNKQpBc20ucDIgPC0gZml4ZWYoTTEpWzFdICsgZml4ZWYoTTEpWzNdICsgcmFuZWYoTTEpWzgsIDFdICMgIkEiIGZvciBNLiBzcGlub3N1bSBpbiB0aGUgcGxvdCAyIChTTSkKQXNtLnAzIDwtIGZpeGVmKE0xKVsxXSArIGZpeGVmKE0xKVszXSArIHJhbmVmKE0xKVsxMSwgMV0gIyAiQSIgZm9yIE0uIHNwaW5vc3VtIGluIHRoZSBwbG90IDMgKFNNKQpBc3MucDEgPC0gZml4ZWYoTTEpWzFdICsgZml4ZWYoTTEpWzRdICsgcmFuZWYoTTEpWzYsIDFdICMgIkEiIGZvciBTLiBmaWxhZ2lub2lkZXMgaW4gdGhlIHBsb3QgMSAoU00pCkFzcy5wMiA8LSBmaXhlZihNMSlbMV0gKyBmaXhlZihNMSlbNF0gKyByYW5lZihNMSlbOSwgMV0gIyAiQSIgZm9yIFMuIGZpbGFnaW5vaWRlcyBpbiB0aGUgcGxvdCAyIChTTSkKQXNzLnAzIDwtIGZpeGVmKE0xKVsxXSArIGZpeGVmKE0xKVs0XSArIHJhbmVmKE0xKVsxMiwgMV0gIyAiQSIgZm9yIFMuIGZpbGFnaW5vaWRlcyBpbiB0aGUgcGxvdCAzIChTTSkKCnBhcihtZnJvdyA9IGMoMiwgMiksIGNleCA9IDAuNzUpICMgKEZpZ3VyZSA0KQojIEdyYXNzZXMgVyBzaXRlCnBsb3QoR1ckdGltZSwgR1ckbGZtYywgeGxhYiA9ICIiLCB5bGFiID0gIkxGTUMgKCUpIiwgeWxpbSA9IGMoMCwgMTAwKSwgbWFpbiA9ICJHcmFzc2VzIFcgc2l0ZSIsIGNvbCA9ICJkYXJrZ3JlZW4iKQpjdXJ2ZSgoKEFndy5wNCAtIFdndykgLyAoMSArIGV4cCgoTSAtIHgpIC8gUykpICsgV2d3KSwgbHdkID0gMSwgbHR5ID0gMiwgYWRkID0gVCwgY29sID0gImRhcmtncmVlbiIpICMgcGxvdCA0CmN1cnZlKCgoQWd3LnA1IC0gV2d3KSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXZ3cpLCBsd2QgPSAxLCBsdHkgPSAyLCBhZGQgPSBULCBjb2wgPSAiZGFya2dyZWVuIikgIyBwbG90IDUKY3VydmUoKChBZ3cucDYgLSBXZ3cpIC8gKDEgKyBleHAoKE0gLSB4KSAvIFMpKSArIFdndyksIGx3ZCA9IDEsIGx0eSA9IDIsIGFkZCA9IFQsIGNvbCA9ICJkYXJrZ3JlZW4iKSAjIHBsb3QgNgojIEdyYXNzZXMgRSBzaXRlCnBsb3QoR0UkdGltZSwgR0UkbGZtYywgeGxhYiA9ICIiLCB5bGFiID0gIkxGTUMgKCUpIiwgeWxpbSA9IGMoMCwgMTAwKSwgbWFpbiA9ICJHcmFzc2VzIEUgc2l0ZSIsIGNvbCA9ICJncmVlbiIpCmN1cnZlKCgoQWdlLnAxIC0gV2dlKSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXZ2UpLCBsd2QgPSAxLCBsdHkgPSAyLCBhZGQgPSBULCBjb2wgPSAiZ3JlZW4iKSAjIHBsb3QgMQpjdXJ2ZSgoKEFnZS5wMiAtIFdnZSkgLyAoMSArIGV4cCgoTSAtIHgpIC8gUykpICsgV2dlKSwgbHdkID0gMSwgbHR5ID0gMiwgYWRkID0gVCwgY29sID0gImdyZWVuIikgIyBwbG90IDIKY3VydmUoKChBZ2UucDMgLSBXZ2UpIC8gKDEgKyBleHAoKE0gLSB4KSAvIFMpKSArIFdnZSksIGx3ZCA9IDEsIGx0eSA9IDIsIGFkZCA9IFQsIGNvbCA9ICJncmVlbiIpICMgcGxvdCAzCiMgTS4gU3Bpbm9zdW0gKEUgc2l0ZSkKcGxvdChTTSR0aW1lLCBTTSRsZm1jLCB4bGFiID0gIlRpbWUgKGRheXMpIiwgeWxhYiA9ICJMRk1DICglKSIsIHlsaW0gPSBjKDAsIDM1MCksIG1haW4gPSAiTS4gc3Bpbm9zdW0iLCBmb250Lm1haW4gPSA0LCBjb2wgPSAiZGFya29yYW5nZSIpCmN1cnZlKCgoQXNtLnAxIC0gV3NtKSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXc20pLCBsd2QgPSAxLCBsdHkgPSAyLCBhZGQgPSBULCBjb2wgPSAiZGFya29yYW5nZSIpICMgcGxvdCAxCmN1cnZlKCgoQXNtLnAyIC0gV3NtKSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXc20pLCBsd2QgPSAxLCBsdHkgPSAyLCBhZGQgPSBULCBjb2wgPSAiZGFya29yYW5nZSIpICMgcGxvdCAyCmN1cnZlKCgoQXNtLnAzIC0gV3NtKSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXc20pLCBsd2QgPSAxLCBsdHkgPSAyLCBhZGQgPSBULCBjb2wgPSAiZGFya29yYW5nZSIpICMgcGxvdCAzCiMgUy4gZmlsYWdpbm9pZGVzIChFIHNpdGUpCnBsb3QoU1MkdGltZSwgU1MkbGZtYywgeGxhYiA9ICJUaW1lIChkYXlzKSIsIHlsYWIgPSAiTEZNQyAoJSkiLCB5bGltID0gYygwLCAzNTApLCBtYWluID0gIlMuIGJyYWN0ZW9sYWN0dXMiLCBmb250Lm1haW4gPSA0LCBjb2wgPSAiZGFya3JlZCIpCmN1cnZlKCgoQXNzLnAxIC0gV3NzKSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXc3MpLCBsd2QgPSAxLCBsdHkgPSAyLCBhZGQgPSBULCBjb2wgPSAiZGFya3JlZCIpICMgcGxvdCAxCmN1cnZlKCgoQXNzLnAyIC0gV3NzKSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXc3MpLCBsd2QgPSAxLCBsdHkgPSAyLCBhZGQgPSBULCBjb2wgPSAiZGFya3JlZCIpICMgcGxvdCAyCmN1cnZlKCgoQXNzLnAzIC0gV3NzKSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXc3MpLCBsd2QgPSAxLCBsdHkgPSAyLCBhZGQgPSBULCBjb2wgPSAiZGFya3JlZCIpICMgcGxvdCAzCmBgYAoKYGBge3J9CiMgMi4yLjMgLSBEZXJpdmF0aXZlcyAoZHJ5aW5nIHNwZWVkKSAtLS0tCmxmbWMuR1cgPC0gZXhwcmVzc2lvbigoQWd3IC0gV2d3KSAvICgxICsgZXhwKChNIC0geCkgLyBTKSkgKyBXZ3cpICMgTEZNQyBkeW5hbWljcyBmb3IgZ3Jhc3NlcyBpbiB0aGUgVyBzaXRlCmRzLkdXIDwtIGRlcml2KGxmbWMuR1csICJ4IiwgZnVuY3Rpb24uYXJnID0gVFJVRSkgIyBEcnlpbmcgc3BlZWQgZm9yIGdyYXNzZXMgaW4gdGhlIFcgc2l0ZQoKbGZtYy5HRSA8LSBleHByZXNzaW9uKChBZ2UgLSBXZ2UpIC8gKDEgKyBleHAoKE0gLSB4KSAvIFMpKSArIFdnZSkgIyBMRk1DIGR5bmFtaWNzIGZvciBncmFzc2VzIGluIHRoZSBFIHNpdGUKZHMuR0UgPC0gZGVyaXYobGZtYy5HRSwgIngiLCBmdW5jdGlvbi5hcmcgPSBUUlVFKSAjIERyeWluZyBzcGVlZCBmb3IgZ3Jhc3NlcyBpbiB0aGUgRSBzaXRlCgpsZm1jLlNNIDwtIGV4cHJlc3Npb24oKEFzbSAtIFdzbSkgLyAoMSArIGV4cCgoTSAtIHgpIC8gUykpICsgV3NtKSAjIExGTUMgZHluYW1pY3MgZm9yIE0uIHNwaW5vc3VtCmRzLlNNIDwtIGRlcml2KGxmbWMuU00sICJ4IiwgZnVuY3Rpb24uYXJnID0gVFJVRSkgIyBEcnlpbmcgc3BlZWQgZm9yIE0uIHNwaW5vc3VtCgpsZm1jLlNTIDwtIGV4cHJlc3Npb24oKEFzcyAtIFdzcykgLyAoMSArIGV4cCgoTSAtIHgpIC8gUykpICsgV3NzKSAjIExGTUMgZHluYW1pY3MgZm9yIFMuIGZpbGFnaW5vaWRlcwpkcy5TUyA8LSBkZXJpdihsZm1jLlNTLCAieCIsIGZ1bmN0aW9uLmFyZyA9IFRSVUUpICMgIyBEcnlpbmcgc3BlZWQgZm9yIFMuIGZpbGFnaW5vaWRlcwoKeHZlYyA8LSBzZXEoMCwgMTAwLCBsZW5ndGggPSAxMDAwKQp5LmRzLkdXIDwtIGRzLkdXKHh2ZWMpCnkuZHMuR0UgPC0gZHMuR0UoeHZlYykKeS5kcy5TTSA8LSBkcy5TTSh4dmVjKQp5LmRzLlNTIDwtIGRzLlNTKHh2ZWMpCmBgYAoKYGBge3J9CnBhcihtZnJvdyA9IGMoMiwgMiksIGNleCA9IDAuNzUpICMgKEZpZ3VyZSA0KQpwbG90KHh2ZWMsIHkuZHMuR1csCiAgeGxpbSA9IGMoMCwgODApLCB5bGltID0gYygwLCA0KSwKICB5bGFiID0gIkRyeWluZyBzcGVlZCIsIHhsYWIgPSAiVGltZSAoZGF5cykiLCB0eXBlID0gIm4iLCBtYWluID0gIkdyYXNzZXMgaW4gdGhlIFcgc2l0ZSIKKQpsaW5lcyh4dmVjLCAtMSAqIChhdHRyKHkuZHMuR1csICJncmFkIikpLCBsd2QgPSAyLCBsdHkgPSAxLCBjb2wgPSAiZGFya2dyZWVuIikKcGxvdCh4dmVjLCB5LmRzLkdFLAogIHhsaW0gPSBjKDAsIDgwKSwgeWxpbSA9IGMoMCwgNCksCiAgeWxhYiA9ICJEcnlpbmcgc3BlZWQiLCB4bGFiID0gIlRpbWUgKGRheXMpIiwgdHlwZSA9ICJuIiwgbWFpbiA9ICJHcmFzc2VzIGluIHRoZSBFIHNpdGUiCikKbGluZXMoeHZlYywgLTEgKiAoYXR0cih5LmRzLkdFLCAiZ3JhZCIpKSwgbHdkID0gMiwgbHR5ID0gMSwgY29sID0gImdyZWVuIikKcGxvdCh4dmVjLCB5LmRzLlNNLAogIHhsaW0gPSBjKDAsIDgwKSwgeWxpbSA9IGMoMCwgNCksCiAgeWxhYiA9ICJEcnlpbmcgc3BlZWQiLCB4bGFiID0gIlRpbWUgKGRheXMpIiwgdHlwZSA9ICJuIiwgbWFpbiA9ICJNLiBzcGlub3N1bSIsIGZvbnQubWFpbiA9IDQKKQpsaW5lcyh4dmVjLCAtMSAqIChhdHRyKHkuZHMuU00sICJncmFkIikpLCBsd2QgPSAyLCBsdHkgPSAxLCBjb2wgPSAiZGFya29yYW5nZSIpCnBsb3QoeHZlYywgeS5kcy5TUywKICB4bGltID0gYygwLCA4MCksIHlsaW0gPSBjKDAsIDQpLAogIHlsYWIgPSAiRHJ5aW5nIHNwZWVkIiwgeGxhYiA9ICJUaW1lIChkYXlzKSIsIHR5cGUgPSAibiIsIG1haW4gPSAiUy4gZmlsYWdpbm9pZGVzIiwgZm9udC5tYWluID0gNAopCmxpbmVzKHh2ZWMsIC0xICogKGF0dHIoeS5kcy5TUywgImdyYWQiKSksIGx3ZCA9IDIsIGx0eSA9IDEsIGNvbCA9ICJkYXJrcmVkIikKYGBgCgpgYGB7cn0KIyAyLjMgLS0tLSBUZXN0aW5nIGRpZmZlcmVuY2VzIGFtb25nIGxlYWYgdHlwZXMgLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLQoKIyAyLjMuMSAtIFBhcmFtZXRlciAiQSIgLS0tLQpjb250cmFzdChlbW1lYW5zKE0xLCB+bGVhZi50eXBlLCBwYXJhbSA9ICJBIiksICJwYWlyd2lzZSIpCmBgYAoKLSAgIE1heGltdW0gTEZNQyBpbiBncmFzc2VzIGZyb20gdGhlIFcgc2l0ZSBpcyBsb3dlciB0aGFuIGdyYXNzZXMgZnJvbSB0aGUgRSBzaXRlIChHcmFzcyBXIC0gR3Jhc3MgRSkKLSAgIE1heGltdW0gTEZNQyBpbiBncmFzc2VzIGlzIGxvd2VyIHRoYW4gc2hydWJzIChHcmFzcyBXIC0gTS4gc3Bpbm9zdW0pLAotICAgTWF4aW11bSBMRk1DIGluIGdyYXNzZXMgaXMgbG93ZXIgdGhhbiBzaHJ1YnMgKEdyYXNzIFcgLSBTLiBicmFjdGVvbGFjdHVzKQotICAgTWF4aW11bSBMRk1DIGluIGdyYXNzZXMgaXMgbG93ZXIgdGhhbiBzaHJ1YnMgKEdyYXNzIEUgLSBNLiBzcGlub3N1bSkKLSAgIE1heGltdW0gTEZNQyBpbiBncmFzc2VzIGlzIGxvd2VyIHRoYW4gc2hydWJzIChHcmFzcyBFIC0gUy4gYnJhY3Rlb2xhY3R1cykKLSAgIFRoZXJlIGFyZSBub3QgZGlmZmVyZW5jZSBiZXR3ZWVuIHRoZSB0d28gc2hydWIgc3BlY2llcyAoTS4gc3Bpbm9zdW0gLSBTLiBicmFjdGVvbGFjdHVzKS4KCmBgYHtyfQojIDIuMy4yIC0gUGFyYW1ldGVyICJ3IiAtLS0tCmNvbnRyYXN0KGVtbWVhbnMoTTEsIH5sZWFmLnR5cGUsIHBhcmFtID0gInciKSwgInBhaXJ3aXNlIikKYGBgCgotICAgTWluaW11bSBMRk1DIGluIGdyYXNzZXMgZnJvbSB0aGUgVyBzaXRlIGlzIGhpZ2hlciB0aGFuIGdyYXNzZXMgZnJvbSB0aGUgRSBzaXRlIChHcmFzcyBXIC0gR3Jhc3MgRSkKLSAgIE1pbmltdW0gTEZNQyBpbiBncmFzc2VzIGlzIGxvd2VyIHRoYW4gc2hydWJzIChHcmFzcyBXIC0gTS4gc3Bpbm9zdW0pCi0gICBNaW5pbXVtIExGTUMgaW4gZ3Jhc3NlcyBpcyBsb3dlciB0aGFuIHNocnVicyAoR3Jhc3MgVyAtIFMuIGJyYWN0ZW9sYWN0dXMpCi0gICBNaW5pbXVtIExGTUMgaW4gZ3Jhc3NlcyBpcyBsb3dlciB0aGFuIHNocnVicyAoR3Jhc3MgRSAtIE0uIHNwaW5vc3VtKQotICAgTWluaW11bSBMRk1DIGluIGdyYXNzZXMgaXMgbG93ZXIgdGhhbiBzaHJ1YnMgKEdyYXNzIEUgLSBTLiBicmFjdGVvbGFjdHVzKQotICAgVGhlcmUgYXJlIG5vdCBkaWZmZXJlbmNlIGJldHdlZW4gdGhlIHR3byBzaHJ1YiBzcGVjaWVzIChNLiBzcGlub3N1bSAtIFMuIGJyYWN0ZW9sYWN0dXMpLgoKIyMgMykgTk9OTElORUFSIEZJWEVELUVGRkVDCgpSZXNwb25zZSBmdW5jdGlvbjogbGZtYyA9IChBIC0gdykgLyAoMSArIGV4cCgobSAtIHRpbWUpL3MpKSkgKyB3CgpgYGB7cn0KIyAzLjEuMSAtIEJhc2Ugbm9ubGluZWFyIG1vZGVsIC0tLS0KbmwuZmUuMCA8LSBubHMobGZtYyB+ICgoQSAtIHcpIC8gKDEgKyBleHAoKG0gLSB0aW1lKSAvIHMpKSkgKyB3LAogIHN0YXJ0ID0gYyhBID0gMTUsIG0gPSAzMCwgcyA9IC0xNywgdyA9IDUpLAogIGRhdGEgPSBsZm1jCikKIyBoZXJlLCB0aW1lIGlzIHRoZSBvbmx5IHByZWRpY3RvciB2YXJpYWJsZQpubC5mZS4wCmBgYAoKYGBge3J9CmIgPC0gY29lZihubC5mZS4wKSAjIGNvZWZmaWNpZW50cyBvZiB0aGUgYmFzZSBtb2RlbApiCmBgYAoKYGBge3J9CiMgMy4xLjIgLSBMZWFmIHR5cGUgZml4ZWQgZWZmZWN0IC0tLS0KbmwuZmUuMSA8LSBubHMobGZtYyB+ICgoQVtsZWFmLnR5cGVdIC0gd1tsZWFmLnR5cGVdKSAvICgxICsgZXhwKChtIC0gdGltZSkgLyBzKSkpICsgd1tsZWFmLnR5cGVdLAogIHN0YXJ0ID0gbGlzdChBID0gcmVwKGJbMV0sIDQpLCBtID0gMjUsIHMgPSAtMTAsIHcgPSByZXAoYls0XSwgNCkpLAogIGRhdGEgPSBsZm1jCikgIyBhbmQgZml0IGlzIG1hZGUgdXNpbmcgdGhlIGJhc2UgbW9kZWwncyBjb2VmZmljaWVudHMgYXMgc3RhcnQgdmFsdWVzCm5sLmZlLjEKYGBgCgpgYGB7cn0KYjEgPC0gY29lZihubC5mZS4xKSAjIGNvZWZmaWNpZW50cyBvZiB0aGUgbW9kZWwgd2l0aCBsZWFmIHR5cGUgYW5kIHRpbWUgYXMgcHJlZGljdG9yIHZhcmlhYmxlcwpiMQpgYGAKCmBgYHtyfQojIDMuMS4zIC0gUGxvdCBmaXhlZCBlZmZlY3QgLS0tLQpubC5mZS4yIDwtIG5scyhsZm1jIH4gKChBW2dyb3VwXSAtIHdbZ3JvdXBdKSAvICgxICsgZXhwKChtIC0gdGltZSkgLyBzKSkpICsgd1tncm91cF0sCiAgc3RhcnQgPSBsaXN0KEEgPSByZXAoYyhiMVsxXSwgYjFbMl0sIGIxWzNdLCBiMVs0XSksIDMpLCBtID0gMjUsIHMgPSAtMTAsIHcgPSByZXAoYyhiMVs3XSwgYjFbOF0sIGIxWzldLCBiMVsxMF0pLCAzKSksCiAgZGF0YSA9IGxmbWMKKSAjIHRoZSBtb2RlbCBpcyBmaXQgdXNpbmcgImIxIiBhcyB0aGUgc3RhcnQgdmFsdWVzCm5sLmZlLjIKYGBgCgpgYGB7cn0KSUNfdGFiKG5sLmZlLjAsIG5sLmZlLjEsIG5sLmZlLjIpCmBgYAoKYGBge3J9CiMgUmVzaWR1YWwgY2hlY2tpbmc6CnBhcihtZnJvdyA9IGMoMSwgMSkpCnBsb3QobmwuZmUuMikgIyBoZXRlcm9jZWRhc3RpY2l0eSBpbiB0aGUgcmVzaWR1YWxzCmBgYAoKYGBge3J9Cmhpc3QocmVzaWQobmwuZmUuMiksIG1haW4gPSAiIiwgeGxhYiA9ICJSZXNpZHVhbHMiKSAjIHNsaWdodGx5IHNrZXcgZGlzdHJpYnV0aW9uCmBgYAoKYGBge3J9CnFxbm9ybShyZXNpZChubC5mZS4yKSkKcXFsaW5lKHJlc2lkKG5sLmZlLjIpKQpgYGAKClRoaXMgbW9kZWwgYXNzdW1lczoKCi0gICAiQSIgYW5kICJ3IiB2YXJ5IHdpdGggbGVhZiB0eXBlIGFuZCBwbG90IChmaXhlZC1lZmZlY3RzKS4KLSAgICJtIiBhbmQgInMiIGFyZSB1bmlxdWUgZm9yIGFsbCB0aGUgcGxvdHMgYW5kIGxlYWYgdHlwZXMuCi0gICBWYXJpYW5jZSBpcyBob21vZ2VuZW91cwotICAgVGhlcmUgaXMgbm8gdGVtcG9yYWwgZGVwZW5kZW5jZQotICAgUmVzaWR1YWxzIGZvbGxvd2luZyBhIG5vcm1hbCBkaXN0cmlidXRpb24KCmBgYHtyfQpNMiA8LSBubC5mZS4yCnN1bW1hcnkoTTIpCmBgYAoKIyMgNCAtIExJTkVBUiBNSVhFRC1FRkZFQ1RTIE1PREVMIChNMykKClJlc3BvbnNlIGZ1bmN0aW9uOiBsZm1jID0gzrIwICsgzrIxeFRpbWUgKyDOsjJ4R0UgKyDOsjN4U00gKyDOsjR4U1MgzrI1eEdFeFRpbWUgKyDOsjZ4U014VGltZSArIM6yN3hTU3hUaW1lCgp3aGVyZSBHRSwgU00sIGFuZCBTUyBhcmUgZHVtbXkgdmFyaWFibGVzICgwIG9yIDEpIGNyZWF0ZWQgdG8gaW5kaWNhdGUgZGlmZmVyZW5jZXMgYmV0d2VlbiB0aGUgYmFzZSBsZXZlbCBvZiAibGVhZiB0eXBlIiwgaS5lLiwgZ3Jhc3NlcyBmcm9tIHRoZSBXIHNpdGUgW0dXXSwgYW5kIHRoZSByZW1haW5pbmcgbGV2ZWxzLgoKIyMjIDQuMS4xIEJhc2UgbGluZWFyIG1peGVkLWVmZmVjdCBtb2RlCgpgYGB7cn0KIyA0LjEuMSAtIEJhc2UgbGluZWFyIG1peGVkLWVmZmVjdCBtb2RlbCAtLS0tCmwubWUgPC0gbG1lKGxmbWMgfiB0aW1lICogbGVhZi50eXBlLCByYW5kb20gPSB+IDEgfCBwbG90LCBkYXRhID0gbGZtYykKCiMgUmVzaWR1YWwgY2hlY2tpbmc6CnBsb3QobC5tZSkgIyBoZXRlcm9jZWRhc3RpY2l0eSBpbiB0aGUgcmVzaWR1YWxzCmBgYAoKYGBge3J9CnBsb3QoQUNGKGwubWUsIHJlc1R5cGUgPSAibm9ybWFsaXplZCIsIG5hLmFjdGlvbiA9IG5hLm9taXQpLCBhbHBoYSA9IDAuMDUsIGdyaWQgPSBUUlVFKSAjIHRlbXBvcmFsIGNvcnJlbGF0aW9uIGlzIG9ic2VydmVkIGF0IGxhZyAxCmBgYAoKYGBge3J9CiMgNC4xLjIgLSBWYXJpYW5jZSBtb2RlbGxpbmcgLS0tLQpsLm1lLnZtIDwtIHVwZGF0ZShsLm1lLCB3ZWlnaHRzID0gdmFySWRlbnQoZm9ybSA9IH4gMSB8IGxlYWYudHlwZSkpCgojIFJlc2lkdWFsIGNoZWNraW5nOgpwbG90KGwubWUudm0pICMgdmFyaWFuY2VzIHRoZXJlIHNlZW0gdG8gYmUgd2VsbCBtb2RlbGVkCmBgYAoKYGBge3J9CkFJQyhsLm1lLCBsLm1lLnZtKSAjIG1vZGVsaW5nIHRoZSB2YXJpYW5jZSBpbXByb3ZlcyB0aGUgbW9kZWwgZml0CmBgYAoKYGBge3J9CiMgNC4xLjMgLSBUZW1wb3JhbCBjb3JyZWxhdGlvbiBtb2RlbGxpbmcgLS0tLQpsLm1lLnZtLnRtIDwtIHVwZGF0ZShsLm1lLnZtLCBjb3JyZWxhdGlvbiA9IGNvckFSTUEocCA9IDEsIHEgPSAxLCBmb3JtID0gfjEpKQpwbG90KEFDRihsLm1lLnZtLnRtLCByZXNUeXBlID0gIm5vcm1hbGl6ZWQiLCBuYS5hY3Rpb24gPSBuYS5vbWl0KSwgYWxwaGEgPSAwLjA1LCBncmlkID0gVFJVRSkgIyB0aGUgdGVtcG9yYWwgZGVwZW5kZW5jZSBpcyBtb2RlbGVkCmBgYAoKYGBge3J9CkFJQyhsLm1lLnZtLCBsLm1lLnZtLnRtKSAjIG1vZGVsaW5nIHRoZSB0ZW1wb3JhbCBjb3JyZWxhdGlvbiBpbXByb3ZlcyB0aGUgbW9kZWwgZml0CmBgYAoKYGBge3J9CnBsb3QobC5tZS52bS50bSkKYGBgCgpgYGB7cn0KcGxvdChBQ0YobC5tZS52bS50bSwgcmVzVHlwZSA9ICJub3JtYWxpemVkIiwgbmEuYWN0aW9uID0gbmEub21pdCksIGFscGhhID0gMC4wNSwgZ3JpZCA9IFRSVUUpCmBgYAoKYGBge3J9Cmhpc3QocmVzaWQobC5tZS52bS50bSwgdHlwZSA9ICJub3JtYWxpemVkIiksIG1haW4gPSAiIiwgeGxhYiA9ICJSZXNpZHVhbHMiKQpgYGAKCmBgYHtyfQpxcW5vcm0ocmVzaWQobC5tZS52bS50bSwgdHlwZSA9ICJub3JtYWxpemVkIikpCnFxbGluZShyZXNpZChsLm1lLnZtLnRtLCB0eXBlID0gIm5vcm1hbGl6ZWQiKSkKYGBgCgpgYGB7cn0KIyA0LjEuNCAtIEV2YWx1YXRpb24gb2YgdGhlIGZpeGVkLWVmZmVjdHMgLS0tLQpsLm1lLm1sIDwtIHVwZGF0ZShsLm1lLnZtLnRtLCBtZXRob2QgPSAiTUwiKSAjIHRoZSBtb2RlbCBpcyBmaXR0ZWQgYnkgTWF4aW11bSBsaWtlbGlob29kCmwubWUubWwuMSA8LSB1cGRhdGUobC5tZS5tbCwgfiAuIC0gdGltZTpsZWFmLnR5cGUpICMgdGhlIG1vZGVsIGlzIHJlZHVjZWQgcmVtb3ZpbmcgdGhlIGludGVyYWN0aW9uICJ0aW1lIHggbGVhZiB0eXBlIgpJQ190YWIobC5tZS5tbCwgbC5tZS5tbC4xKSAjIHRoZSBBSUMgY29tcGFyaXNvbiBzdWdnZXN0cyB0aGUgaW50ZXJhY3Rpb24gdGVybSB0byBiZSBpbXBvcnRhbnQKYGBgCgpgYGB7ciwgcGFnZWQucHJpbnQ9RkFMU0V9Ck0zIDwtIGwubWUudm0udG0Kc3VtbWFyeShNMykKYGBgCgpgYGB7cn0KTTQgPC0gbG0obGZtYyB+IHRpbWUgKiBncm91cCwgZGF0YSA9IGxmbWMpCnN1bW1hcnkoTTQpCmBgYAoKYGBge3J9CiMgUmVzaWR1YWwgY2hlY2tpbmc6CnBhcihtZnJvdyA9IGMoMiwgMikpCnBsb3QoTTQpICMgYSBjbGVhciBwYXR0ZXJuIGluIHRoZSByZXNpZHVhbHMgaXMgb2JzZXJ2ZWQgKG1vZGVsIGFzc3VtcHRpb25zIGFyZSBub3QgbWV0KQpgYGAKCmBgYHtyfQojIDYgLSBOVUxMIE1PREVMIChNNSkgIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMKCk01IDwtIGxtKGxmbWMgfiAxLCBkYXRhID0gbGZtYykKc3VtbWFyeShNNSkKYGBgCgpgYGB7cn0KSUNfdGFiKE0yLCBNMywgTTQsIE01KSAjICwgd2VpZ2h0cyA9IFQsIGRlbHRhID0gVFJVRSwgYmFzZT1ULCBzb3J0ID0gVFJVRSkgIyB0aGUgbm9ubGluZWFyIG1peGVkLWVmZmVjdHMgbW9kZWwgaXMgY2xlYXJseSB0aGUgYmVzdApgYGAKCmBgYHtyfQpJQ190YWIoTTEpCmBgYAoKYGBge3J9CnBhcihtZnJvdyA9IGMoMiwgMiksIGNleCA9IDAuNjUpICMgKEZpZ3VyZSA1KQp2MSA8LSBjKDEsIDIsIDMsIDQpCnYyIDwtIGMoIkdXIiwgIkdFIiwgIlNNIiwgIlNTIikKCiMgTm9ubGluZWFyIG1peGVkLWVmZmVjdHMgbW9kZWwgKE0xKQpyTTEgPC0gcmVzaWQoTTEsIHR5cGUgPSAibm9ybWFsaXplZCIpCmZNMSA8LSBmaXR0ZWQoTTEpCnBsb3QoZk0xLCByTTEsIHhsYWIgPSAiRml0dGVkIExGTUMgKCUpIiwgeWxhYiA9ICJSZXNpZHVhbHMiLCB5bGltID0gYygtMywgMykpCmFibGluZShhID0gMCwgYiA9IDAsIGx3ZCA9IDIpCmJveHBsb3Qock0xIH4gbGZtYyRsZWFmLnR5cGUsIHlsYWIgPSAiUmVzaWR1YWxzIiwgeGxhYiA9ICJMZWFmIHR5cGUiLCB4YXh0ID0gIm4iLCB5bGltID0gYygtMywgMykpCmF4aXMoc2lkZSA9IDEsIGF0ID0gdjEsIGxhYmVscyA9IHYyLCBsYXMgPSAxKQpwbG90KHJNMSB+IGxmbWMkdGltZSwgeWxhYiA9ICJSZXNpZHVhbHMiLCB4bGFiID0gIlRpbWUgKGRheXMpIiwgeWxpbSA9IGMoLTMsIDMpKQphYmxpbmUoYSA9IDAsIGIgPSAwLCBsdyA9IDIpCnFxbm9ybShyTTEsIG1haW4gPSAiIiwgeGxhYiA9ICJUaGVvcmV0aWNhbCBxdWFudGlsZXMiLCB5bGFiID0gIlNhbXBsZSBxdWFudGlsZXMiKQpxcWxpbmUock0xKQpgYGAKCmBgYHtyfQojIE5vbmxpbmVhciBmaXhlZC1lZmZlY3RzIG1vZGVsIChNMikKcGFyKG1mcm93ID0gYygyLCAyKSwgY2V4ID0gMC42NSkgIyAoRmlndXJlIDUpCgpyTTIgPC0gcmVzaWQoTTIpCmZNMiA8LSBmaXR0ZWQoTTIpCnBsb3QoZk0yLCByTTIsIHhsYWIgPSAiRml0dGVkIExGTUMgKCUpIiwgeWxhYiA9ICJSZXNpZHVhbHMiKQphYmxpbmUoYSA9IDAsIGIgPSAwLCBsd2QgPSAyKQpib3hwbG90KHJNMiB+IGxmbWMkbGVhZi50eXBlLCB5bGFiID0gIlJlc2lkdWFscyIsIHhsYWIgPSAiTGVhZiB0eXBlIiwgeGF4dCA9ICJuIikKYXhpcyhzaWRlID0gMSwgYXQgPSB2MSwgbGFiZWxzID0gdjIsIGxhcyA9IDEpCnBsb3Qock0yIH4gbGZtYyR0aW1lLCB5bGFiID0gIlJlc2lkdWFscyIsIHhsYWIgPSAiVGltZSAoZGF5cykiKQphYmxpbmUoYSA9IDAsIGIgPSAwLCBsdyA9IDIpCnFxbm9ybShyTTIsIG1haW4gPSAiIiwgeGxhYiA9ICJUaGVvcmV0aWNhbCBxdWFudGlsZXMiLCB5bGFiID0gIlNhbXBsZSBxdWFudGlsZXMiKQpxcWxpbmUock0yKQpgYGAKCmBgYHtyfQpwYXIobWZyb3cgPSBjKDIsIDIpLCBjZXggPSAwLjY1KQoKIyBMaW5lYXIgbWl4ZWQtZWZmZWN0cyBtb2RlbCAoTTMpCnJNMyA8LSByZXNpZChNMywgdHlwZSA9ICJub3JtYWxpemVkIikKZk0zIDwtIGZpdHRlZChNMykKcGxvdChmTTMsIHJNMywgeGxhYiA9ICJGaXR0ZWQgTEZNQyAoJSkiLCB5bGFiID0gIlJlc2lkdWFscyIsIHlsaW0gPSBjKC0zLCAzKSkKYWJsaW5lKGEgPSAwLCBiID0gMCwgbHdkID0gMikKYm94cGxvdChyTTMgfiBsZm1jJGxlYWYudHlwZSwgeWxhYiA9ICJSZXNpZHVhbHMiLCB4bGFiID0gIkxlYWYgdHlwZSIsIHhheHQgPSAibiIsIHlsaW0gPSBjKC0zLCAzKSkKYXhpcyhzaWRlID0gMSwgYXQgPSB2MSwgbGFiZWxzID0gdjIsIGxhcyA9IDEpCnBsb3Qock0zIH4gbGZtYyR0aW1lLCB5bGFiID0gIlJlc2lkdWFscyIsIHhsYWIgPSAiVGltZSAoZGF5cykiLCB5bGltID0gYygtMywgMykpCmFibGluZShhID0gMCwgYiA9IDAsIGx3ID0gMikKcXFub3JtKHJNMywgbWFpbiA9ICIiLCB5bGFiID0gIlNhbXBsZSBxdWFudGlsZXMiLCB4bGFiID0gIlRoZW9yZXRpY2FsIHF1YW50aWxlcyIpCnFxbGluZShyTTMpCmBgYAoKYGBge3J9CnBhcihtZnJvdyA9IGMoMiwgMiksIGNleCA9IDAuNjUpCgojIENsYXNpY2FsIHJlZ3Jlc3Npb24gKE00KQpyTTQgPC0gcmVzaWQoTTQpCmZNNCA8LSBmaXR0ZWQoTTQpCnBsb3QoZk00LCByTTQsIHhsYWIgPSAiRml0dGVkIExGTUMgKCUpIiwgeWxhYiA9ICJSZXNpZHVhbHMiKQphYmxpbmUoYSA9IDAsIGIgPSAwLCBsd2QgPSAyKQpib3hwbG90KHJNNCB+IGxmbWMkbGVhZi50eXBlLCB5bGFiID0gIlJlc2lkdWFscyIsIHhsYWIgPSAiIiwgeGF4dCA9ICJuIikKYXhpcyhzaWRlID0gMSwgYXQgPSB2MSwgbGFiZWxzID0gdjIsIGxhcyA9IDEpCnBsb3Qock00IH4gbGZtYyR0aW1lLCB5bGFiID0gIlJlc2lkdWFscyIsIHhsYWIgPSAiVGltZSAoZGF5cykiKQphYmxpbmUoYSA9IDAsIGIgPSAwLCBsdyA9IDIpCnFxbm9ybShyTTQsIG1haW4gPSAiIiwgeGxhYiA9ICJUaGVvcmV0aWNhbCBxdWFudGlsZXMiLCB5bGFiID0gIlNhbXBsZSBxdWFudGlsZXMiKQpxcWxpbmUock00KQpgYGAKCiMgLS0tLQo=